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Abstract 

Detecting New Planets in Transiting Systems 

Jason Steffen 

Chair of the Supervisory Committee: 
Professor Eric Agol 
Physics 

I present an initial investigation into a new planet detection technique that uses 
the transit timing of a known, transiting planet. The transits of a solitary planet 
orbiting a star occur at equally spaced intervals in time. If a second planet is present, 
dynamical interactions within the system will cause the time interval between transits 
to vary. These transit time variations can be used to infer the orbital elements of the 
unseen, perturbing planet. I show analytic expressions for the amplitude of the transit 
time variations in several limiting cases. Under certain conditions the transit time 
variations can be comparable to the period of the transiting planet. 1 also present the 
application of this planet detection technique to existing transit observations of the 
TrES-1 and HD209458 systems. While no convincing evidence for a second planet 
in either system was found from those data, I constrain the mass that a perturbing 
planet could have as a function of the semi- major axis ratio of the two planets and 
the eccentricity of the perturbing planet. Near low-order, mean-motion resonances 
(within about 1% fractional deviation), I find that a secondary planet must generally 
have a mass comparable to or less than the mass of the Earth-showing that these 
data are the first to have sensitivity to sub Earth-mass planets orbiting main sequence 
stars. These results show that TTV will be an important tool in the detection and 
characterization of extrasolar planetary systems. 
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Chapter 1 
INTRODUCTION 



Since the first discovery, a decade ago, of a planet orbiting a distant main-sequence 
star Q]; a new field of astronomy has emerged with the potential to address funda- 
mental questions about our own solar system since we can now compare it with 
other planetary systems. Several search techniques have been employed to iden- 
tify additional planets ( "exoplanets" ) in orbit around distant stars. The primary 
mode for discovery of exoplanets has been the measurement of the stellar radial ve- 
locities via the Doppler effect upon the spectral lines of the host star. Currently 
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planetary microlensing 
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and direct imaging 



171 l24l lid . |27[ . Recently, a large number of planetary transit searches are bein 



carried out which are starting to yield an handful of giant planets |l8, 
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and many more planned searches should reap a large harvest of transiting planets in 
the near future 30|. Despite these successes, the discovery of "terrestrial" extrasolar 
planets, similar in size and mass to the Earth, awaits further developments in all of 
these planet detection techniques. 

In this dissertation I present much of the founding research of a recent addition 
to the repertoire of planet search techniques which consists of looking for additional 
planets in a known, transiting system by analyzing the variation in the time between 
planetary transits. These transit timing variations (TTV), which are caused by the 
gravitational perturbation of a secondary planet, can be used to constrain the orbital 
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elements of the unseen, perturbing planet — even if its mass is comparable to the mass 
of the Earth 3, I, 2^- For the near term, this new technique can probe for planets 
around main-sequence stars that are too small to detect by any other means. This 
arly manifest near mean-motion resonances, which recent works 
and [221 suggest might be very common. 



sensitivitv is particu 



by 



The application of TTV depends upon the discovery and monitoring of transiting 



planetary systems. The first successful detection of a planetary transit [18|, |28[ was for 
a planet that had originally been identified from the reflex motion of the primary star, 
HD209458b 3- 

Since the mass of the planet is degenerate with orbital inclination, 
the planetary status of the companion was confirmed once the transits were seen 
(since planetary transits imply that the system is seen edge-on). i/S" T observations 
of the HD209458 yielded precision measurements of the transit lightcurve [isj] and 
made this the surest planetary candidate around a main sequence star (other than 
our own). 

The first extrasolar planet to be discovered primarily from transit data was the 
OGLE-TR-56b system reported by js^. Existing and future searches for planetary 
transits, such as the COROT jl^, XO j47|, and Kepler missions, are expected to 
identify many, possibly hundreds, of transiting systems. Each of the systems that will 
be discovered is a candidate for an analysis of the timing of the planetary transits. 
These analyses may yield additional planetary discoveries or may provide important 
constraints on the presence of additional bodies in each system. One such system is 
the TrES-1 planetary system, reported by j4|, which was also discovered via planetary 
transits. In a recent paper by David Charbonneau and collaborators the detection 
of thermal emission from the surface of TrES-lb was announced. This paper includes 
the timing of 12 planetary transits. In chapter 5, I present the results from the first 
published TTV analysis of transit data that was designed to identify and constrain 
potential companion planets in such a system. The analysis of these data proved to 



be the first probe for sub earth-mass planets in orbit about a distant, main-sequence 
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star. 

Beyond planet detection, TTV has other apphcations that are important for ex- 
trasolar planet research. For example, while the ratio of the planetary radius to the 

r~i 

stellar radius of a transiting system can be measured with extreme precision '^flj] , the 
absolute radii remain uncertain due to a degeneracy that exists between the radius 
and the mass of the host star — an increase in the mass and radius of the star can 
yield an identical lightcurve and period. This mass-radius degeneracy may be broken 
via TTV provided there is an additional planet in the system. About 10 per cent of 
the stars with known planetary companions have more than one planet, while possibly 
as much as 50 per cent of them show a trend in radial velocity indicative of additional 
planets 2^- If both of the planets is transiting, dynamical interactions be- 



18, m. A measurement of 



tween the planets will alter the timing of the transits |22 , 
these timing variations, combined with radial velocity data, can break the mass-radius 
degeneracy. 

Other applications of TTV have significant ramifications for various theories of 
the formation and evolution of multiple planet systems. One example is that TTV 
may be used to identify the relative inclinations of planetary systems. A sample of 
relative inclinations of planetary systems would constrain the mechanisms by which 

n 

the eccentricities of planets can grow 58]. A second example is that TTV is well 
suited to detect the presence of small, rocky planets that may be trapped in low- 
order, mean-motion resonances. The presence of these close-in, resonant, terrestrial 
planets favors a sequential-accretion model of planet formation over a gravitational 
instability model 72j |. 

Given the multiple motivations of detecting terrestrial planets, breaking the mass- 
radius degeneracy, and constraining planetary system formation and evolution theo- 
ries, the results of this work, along with those of planned developments, should prove 
to be an important tool for this new field of extrasolar planetary science. In this 
dissertation I present analytic and numerical results for transit timing variations due 
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to the presence of a second planet in chapter 2. That chapter shows the transit timing 
variations that are expected in several limiting cases such as non-interacting planets, 
an eccentric exterior perturbing planet with a large period relative to the transiting 
planet, the general transit timing differences for two planets with circular, co-planar 
orbits, the case of exact mean-motion resonance, and the case of two eccentric planets. 
I also present the results from numerical simulations of several known multi-planet 
systems (though they are not transiting). 

In chapter 3, I discuss in more detail some of the applications of the TTV tech- 
nique. I present the results of a preliminary study to characterize the efficiency with 
which the TTV technique can be applied to discover secondary planets in transiting 
systems in chapter 4. Chapters 5 and 6 contain results of TTV analyses of the transit 
times of the known, transiting systems TrES-1 and HD209458. Finally, I give some 
concluding remarks in chapter 7. 
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Chapter 2 



THE TTV SIGNAL 



In this chapter I present various mechanisms that can cause variations in the 
period of a transiting planet. In general, all of these mechanisms are present in a 
physical system. However, a study of several limiting cases serves both to identify 
important relationships among the orbital elements in relevant systems and to identify 
the effects of each individual means of perturbation. Much of this discussion can be 
found in 

I will begin by outlining the coordinate system that I use for these analytic treate- 
ments. I then discuss, in turn, the effects of: noninteracting planets f ^2.2j) . a perturb- 
ing planet on a wide, eccentric orbit f ^2.3p . non- resonant planets with perturbatively 
small eccentricities f ^2.4|) . resonant systems f ^2.5|) . and orbits with arbitrary eccen- 
tricities f ^2.6p . I note that this treatment addresses only systems where the orbits 
are coplanar. 

Throughout the rest of this work I characterize the strength of transit timing 
variations as follows. For a series of transit times, tj, 1 fit the times assuming a 
constant period, P and compute the standard deviation, a, of the difference between 
the nominal and actual times. Mathematically, 




1 1/2 



(2.1) 



where P and to are chosen to minimize a. If the variations are strictly sinusoidal, 
then the amplitude of the timing deviation is simply ^/2 times larger than a. 
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2.1 The Coordinate Systems 

For this study of three body systems the positions of the objects, with arbitrary origin, 
are given by the Cartesian coordinates Rj, i = 0,1,2 where the labels correspond 
to the central body and the two planetary companions respectively. According to 
Newton's law of gravity, the accelerations of the masses are given by 

^-=E^^^ i^r^v (2-2) 



By multiplying the equation for each particle by that particle's mass and adding 
them together, one finds: 

moRo + miRi + m2R2 = 0. (2.3) 

which states that the center of mass of the system experiences no external forces. 
Since light travel time and parallax effects are negligible (see ^3.3. HI . the timing of 
each transit is unaffected by the total velocity or position of the center of mass. So, 
setting 

Rem = = (2.4) 

reduces the differential equations of motion to two, which I take to be that of the two 
planets, Ri and i?2- 

When I solve these equations of motion numerically I use this Cartesian coordi- 
nate system. However, for an analytic investigation it is more convenient to write 



the problem in Jacobian coordinates, a coordinate 



49 



system that is commonly used in 



38, 



For the three-body prob- 



perturbation theory for many bodies [see, e.g. 
lem, the Jacobian coordinates amount to three new coordinates which describe (a) the 
center of mass of the system; (b) the relative position of inner planet and the star (the 
"inner binary"); and (c) the relative position of the outer planet and the barycenter of 
the inner binary (the "outer binary"). To distinguish from the Cartesian coordinates, 
I denote the Jacobian coordinates with a lower case rj. The Jacobian coordinates are 

To = Rem = 0, 
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ri — Ri — Ro, 

r2 = R2 [moRo + miRi]. (2.5) 

mo + mi 

Using fii = mi/M ~ mi/m^, where M = Yl'i=o'^i^ equations of motion ()2.2|) may 

be rewritten in Jacobian coordinates, 

Gmo ri ri - rsi r2i 

ri = T - GMfi2i R - GMfi2^, 

1-fXirf |ri-r2i|^ r^i 

Gmo r2i r2i — ri 

rs = :t--GMij,i- -, (2.6) 

1-At2'^2i |r2i-ri|'^ 

where r2i = /iiri + r2 = R2 - Ro- 

2.2 An interior, noninteracting planet 

The first, hmiting case that I study assumes that the interaction between the planets 
is neghgible. At first glance it may seem that such an investigation would be fruitless 
and, indeed, if the interior planet is the one that transits there are no variations in the 
transit times due to a noninteracting secondary planet. However, the motion of the 
inner binary causes the position of the star to oscillate with respect to the stationary 
center of mass of the system. This oscillation will, in turn, cause the transits of the 
outer planet to vary as illustrated in figure 12.11 

To quantify these deviations, we take the limit as 1^1,^2 in equation (j2.6|) 
which gives the equations of motion 

ri = -Gmo^, 
rl 

V2 = -Gmo^. (2.7) 

^2 

This approximation requires that the periapse of the outer planet be much larger 
than the apoapse of the inner planet, (1 — 62)02 ^ (1 + ei)ai where 01,02 are the 
semi-major axes of the inner and outer binary and 61,62 are the eccentricities. The 
solution to these equations are standard Keplerian trajectories. 

The simplest case to consider is that in which both the inner and outer binary are 
on approximately circular orbits. The transit occurs when the outer planet is nearly 
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9 



5 { 



early 



late 



on time 



Figure 2.1: Cartoon showing how an interior planet that does not directly interact 
with an exterior, transiting planet can change the timing of transits of the exterior 
planet. 
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aligned with the barycenter of the inner binary and its motion during the transit is 
essentially transverse to the hne of sight. The inner planet displaces the star from 
the barycenter of the inner binary by an amount 

xo = -aini sin [27r(i - to) /Pi] (2.8) 

where the inner binary undergoes a transit at time to and Pi is the orbital period of 
the inner binary. Thus, the timing deviation of the mth transit of the outer planet is 

St, ^ ^ ^ _ P2aiHism[2nimP,-to)/Pi] 

V2 - vo 27ra2 

where Vi is the velocity of the ith body with respect to the line of sight. Typically 
Vo <^ V2, so we have neglected vo in the second expression in the previous equation. 
The standard deviation of timing variations over many orbits is 

- = «^'^«"' = S^- 

Note that if the periods have a ratio P2:Pi of the form q:l for some integer q, then 
the perturbations disappear because the argument of the sine function is the same 
each orbit. 

More interesting variations occur if either or both planets are on eccentric orbits. 
Because both planets are following Keplerian orbits, the transit timing variations 
and duration variations can be computed by solving the Kepler problem for each 
Jacobian coordinate. Since wc arc assuming that the planets arc coplanar and edge- 
on, 4 coordinates for each planet suffice to determine the planetary positions: ei,2, 
fli,2, ^1,2 (longitude of pericenter), and /i_2 (true anomaly). As in the circular case, 
the change in the transit timing is approximately St, ~ Xq/v,- 

The position of the star with respect to the barycenter of the inner binary is 



Xq = -Aiirisin[/i + G7i]. 



(2.11) 
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If ai <^ a2, the outer planet is in nearly the same position at the time of each transit 

and its velocity perpendicular to the line of sight is 

27ra2 (1 + 62 cos CiJo) 

V2 = — — (2.12) 

P9 



where we have used the fact that /2 = —tU2 at the timing of the transit. Thus, to 
first order in ai/a2 



g^^^_^2^^vl^^^HJl^^liV---2_ (2.13) 



P2f iiri sin [fi + - ej 

27102(1 + 62 cos W2] 

The standard deviation of 6t2 can be found analytically as well. To do this, I first 
calculate the mean transit deviation (612). Over many transits by the outer planet, 
the inner binary's position populates all of its orbit provided the planets do not have 
a period ratio that is the ratio of two integers. Consequently, I find the mean transit 
deviation by averaging over the probability that the inner binary is at any position 
in its orbit, p{fi) = "n-i/Zi (where rii = 27r/Pi), times the transit deviation at that 
point. This gives 



3 cti . ^ A\ 

= — -/ii — CismciJi. (^•14) 
2 f 2 

The linear scaling with ci is because the star spends more time near apoapse which 
displaces the average position of the star, and hence the transit geometry, away from 
the barycenter of the system. The symmetry of the orbit about w = and vr explains 
the dependence on sintui. 

A similar calculation gives (^^2) resulting standard deviation is 

a2 = {{6tl) - {6t2rY^' 



^^"^^^^^ ri-|(l + cos^t.0l'''. (2.15) 



2^/^7ra2(l + 62 cos W2) 
This agrees with equation ()2.10|) in the limit ei ^ 0. Averaging again over wi and 
tU2, gives 

_P2ai/ii[l-|e?]'/' 
(^2)^,,^. - 23/%a2(l-ei)V4 • ^2-^^) 
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Note that an eccentric inner orbit reduces ct2 because the inner binary spends more 
time near apoapse as the eccentricity increases which, in turn reduces the variation in 
the position of the star when averaged over time. As ei approaches unity for an orbit 
viewed along the major axis, 02 reduces to zero because the minor axis approaches 
zero and there is no resulting variation in the xq position. 

2.3 An exterior planet on a large eccentric orbit 

In this section I include planet-planet interactions to find the timing variations of an 
interior transiting planet — on an orbit with negligible eccentricity — that are caused by 
a perturbing planet that is on an eccentric orbit with a semi-major axis that is much 
larger than that of a transiting planet. In this limit, resonances are not important and 
the ratio of the semi-major axes can be used as a small parameter for a perturbation 
expansion. A general formula for this case has been derived by ^. Here I present a 
shorter derivation which clarifies the primary physical effects for coplanar orbits. 

The essential physics of this scenario is that the presence of the outer planet alters 
the effective mass, and thus the period, of the inner binary — the outer planet acts, in 
the lowest multipole, as a sphere surrounding the inner binary. As the outer planet 
moves along its eccentric trajectory, its proximity to the inner binary changes. The 
periodically changing distance between the perturbing planet and the inner binary 
causes the effective mass of the inner binary to change in tandem with the position 
of the outer planet. So, as the outer planet moves inward, the inner binary slows in 
its orbit; as the outer planet moves outward, the inner binary orbits more rapidly. 
These variations in the period of the inner binary translate to the timing deviations 
of the interior transiting planet. 

The naive model of a pulsating sphere provides the proper scaling relations for the 
dispersion in the timing deviations. I first derive this quantity using this approach 
since it encapsulates the relevant physics of the system. Later, I provide a more 
detailed derivation of the transit timing variations that can be applied more generally. 
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Kepler's Laws state that the period of the inner binary scales as Pi ~ Xj \/Gp\ 
where p\ = mQ/a\ is roughly the average density of material inside the orbit of the 
transiting planet. The outer planet, located a distance r2 from the barycenter, acts 
to reduce this density, now denoted p, and alter the period of the inner binary. The 
difference between the nominal period and the period P that includes the effects of 
the outer planet is 



This result gives the timing deviation that results from one orbit of the transiting 
planet. These deviations add over several orbits — roughly P2/-P1 orbits — which serves 
to increase the typical timing deviation to an amount that is larger than that of a 
single orbit. If we characterize the distance to the second planet as r2 ~ (1 + e)a2 the 
resulting dispersion of timing deviations is 



1 note that this result does not change if the nominal period Pi in the derivation 
includes the average effect of the perturbing planet. 

The item of fundamental importance is the variation in the distance r2 in (j2.17p 
which result in the factor of 62 in equation ()2.18|) . The factor of P2 encapsulates both 
the period of the transiting planet and the time scale over which the timing deviations 
accumulate, the mass ratio gives the relative strengths of the primary, central force 
which causes the orbit of the transiting planet and the secondary, perturbing force 
which causes the deviations in the orbit, and the factor of e2(ai/a2)^ characterizes 
the change in the volume of the orbit of the perturbing planet (volume being appro- 
priate since the lowest multipole approximation is to treat the perturbing planet as 




(2.17) 




3 



(2.18) 
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+ Oir^jr^f. (2.20) 



a spherically symmetric change in the density of the system). 

I now derive, in more detail, the individual transit timing variations as well as the 
dispersion in those deviations that occur in this scenario. The equations describing 
the inner binary can be divided into a Keplerian equation ()2.7|) and a perturbing force 
proportional to m^. The perturbing acceleration bi\ on the inner binary due to the 
outer planet is given by 

5ri = -GM/i2-^^^^^ - GM/is^. (2.19) 
|ri-r2i|^ r^i 

We expand this in a Legendre series and keep terms up to first order in the ratio of 
the radii of the inner and outer orbit, 

e ■■ r „ri ■ r2 

dri = ^ ri-3 — ^ 

ri L r^ 

To compute the perturbed orbital period we must find the change in the force on 
the inner binary due to the outer planet averaged over the orbital period of the inner 
binary. Since the inner binary is nearly circular, the angle of the inner binary is given 
by = /i + roi ^ ni{t — Ti)wi, where I have approximated ei ~ 0. Differentiating 
this with respect to time gives 

0^ = ni{t - Ti) + rii - Tiifi. (2.21) 

Now, we write hi = — 3ni/(2ai)ai, and express cli, Wi, and fi in terms of the radial, 
tangential, and normal components of the force [see section 2.9 of 3|. Plugging these 
expressions into 6i gives a cancellation of most terms to lowest order in ei, and after 
setting the normal force to zero leaves the remaining term 

Oi =ni(l- ^. ^ r) , (2.22) 

V G(mo + mi) J 

where R is the radial disturbing force per unit mass, R = (5ri • ri) = |Gm2ai/r2. 

Thus, the change in the effective mass of the inner binary due to the presence of 
the second planet is actually — |m2(ai/r2)'^, which results in a slight increase in the 
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period of the orbit. The increase in period would be constant if the second planet 
were on a circular orbit. However, for an eccentric orbit, the time variation of r2 
induces a periodic change in the orbital frequency of the inner binary with period 
equal to P2. 

Now, the time of the (A^ + l)th transit occurs at 

»/0+27rAf 



fo 



/o+27rAf 



dfiTi^' 



fo 



1 + 



m2 



mo + mi \r2 



ai 



(2.23) 



where fo is the true anomaly of the inner binary at the time of the first transit. 



Following 



the variable of integration can be changed from /i to /2, the true 



anomaly of the outer planet. 



dfi 



;df2. 



[2.2A) 



Since we are assuming that the orbit of the outer planet is eccentric, r2 = 02(1 — 
62)7(1 + 62 cos /2) , which gives the transit time 

m2 (1 - Co 



tpr. — tn 



NP^ + 



(/2 + 62 sin/2 



(2.25) 



27r(mo + mi) P2 

where /2 is the true anomaly of the outer binary at the timing of the (A^ — 1) transit. 
The unperturbed /2 includes the mean motion, n2(t — T2), which grows linearly with 
time. To find the deviation of the time of transits from a uniform period, I subtract 
off this mean motion as well as NPi which results in 



5ti 



/3 (1 - el) 



-3/2 



m2 



[/2 -n2{t -T2) + 62 sin /2] 



(2.26) 



27r(mo + mi) P2 ' 

This agrees with the expression of jsf in the limit 1 = 1 (i.e. coplanar orbits). 

Numerical calculation of the 3-body problem show that this approximation works 
extremely well in the limit ri -C r2 (see Figure If Pi ^ P2 and the period ratio is 
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Figure 2.2: Contour plot of the logarithm of the dispersion of the normalized timing 
variations, log (aini/ig for an inner circular planet and an outer eccentric planet 
(for example, at the -2.00 contour an orbit lasting 27r years with a perturbing planet 
of mass 1O~^M0 would have a transit timing standard deviation of 10~^ years, or 5 
minutes). The dotted line is the approximation given in equation ()2.28|) . 
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non-rational, then over a long time the transits of the inner planet sample the entire 
phase of the outer planet. Thus, we can find the standard deviation of the transit 
timing variations as in equation ()2.15p 



which is accurate to better than 2 per cent for all 62 and agrees with equation (j2.18|) 
in the limit as 62 — > 0. Figure EIH shows a comparison of this approximation with the 
exact numerical results averaged over W2 (since there is a slight dependence on the 
value of W2). This approximation breaks down for 02(1 — 62) < 5ai since resonances 
and higher order terms contribute strongly when the planets have a close approach. It 
also breaks down for 62 ^ 0.05 since the perturbations to the semi-major axes caused 
by interactions of the planets contribute more strongly than the tidal terms which 
become weaker with smaller eccentricity. 

2.4 Non resonant planets on initially circular orbits 

In this section I calculate the amplitude of timing variations for two planets on nearly 
circular orbits. That is, on orbits where a first-order, perturbative expansion of the 
planetary motion is valid. The resonant forcing terms are most important factors 
that determine the amplitude of the periodic timing variations, even for non-resonant 
planets. The planets interact most strongly at conjunction, so the perturbing planet 
causes a radial kick to the transiting planet inducing eccentricity into its orbit. This 
change in eccentricity, in turn, causes a change in the semi-major axis by an amount 
given by the Tisserand relation 




(2.27) 



since {6ti) = 0, where p{f2) = ^2/ f2- This integral is intractable analytically, but an 
expansion in 62 yields 





da _ 2a5/2g 
~ a3/2 - 1 



(2.29) 
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where a is, in this case, the ratio of the semi-major axes of the planets. 

Both the change in eccentricity and the change in semi-major axis contribute to 
the variations in the period of each planet. This can be seen by Taylor-expanding the 
time-dependent angular velocity of the planet 9 to first order in e 

Tl fl ~1~ 6 COS /")^ 

9 = ^ no + 5n + 2eno cos [A - ro], (2.30) 

(1 — e^) ' 

where no is the nominal mean-motion of the planet, 5n reflects the variations caused 
by changes to the mean-motion, and the oscillatory term of amplitude proportional 
to e quantifies the variations caused by changes to the eccentricity. I first derive the 
timing deviations for when the change in mean-motion dominates the transit time 
variations. Then, I derive an approximation to the variations that result when the 
induced eccentricity dominates the timing deviations — which happens to occur farther 
from the resonances than the mean-motion dominated variations. Finally, I derive 
in detail the eccentricity dominated variations using perturbation theory. I show the 
variations that occur when the planets are in mean-motion resonance in the next 
section. 

2.4-1 Mean-motion dominated variations 

As previously stated, the radial perturbing force on the transiting planet is largest 
when the planets are in conjunction. Since, in this limit where changes to the mean 
motion dominates the timing variations, the planets are not exactly on resonance, 
the longitude of conjunction will drift with time. Eventually, the effects of the radial 
kicks cancel after the longitude drifts by ~ tt in the inertial frame. Thus, the total 
amplitude of the eccentricity grows over a time equal to half of the period of circulation 
of the longitude of conjunction. The closer the planets are to a resonance, the longer 
the period of circulation and thus the larger the eccentricity grows. The change in 
eccentricity causes a change in the semi-major axis and mean motion — the effect we 
study here. Identifying the regime where the timing deviations are dominated by the 
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mean motion will be easiest to identify once both the mean-motion effects and the 
eccentricity effects have been independently quantified and I defer to the end of the 
section B.4.21 for the location of this transition. 

For two planets that are on circular orbits near a j:j + 1 resonance, conjunctions 
occur every Pconj = 27r/(ni — ^2) — jP (I take the limit of large j and we ignore 
factors of order unity). Let us define the fractional distance from resonance, e = 
|1 — (1 + j^^)Pi/ < 1, where e = indicates exact resonance. The longitude of 
conjunction changes with each successive conjunction and ultimately returns to its 
initial value over a period Pcyc = Pj^^^^^- The number of conjunctions per cycle is 
Nc = Pcyc/ Pconj — j "^^^^- Each coujuuction changes the eccentricity of the planets by 
Ae ~ fipertia / (using the perturbation equations for eccentricity and the impulse 
approximation, where fipert is the planet-star mass ratio of the perturbing planet). 
Over half a cycle the eccentricities grow to about NcAe ~ /Xperi(l — Pi/ P2)^^{j^)^^ — 

To find the change in the transit timing that is caused by the 5n term in equation 
I2.3UI I apply the Tisserand relation to the lighter planet (now using subscripts "light" 
and "heavy"), resulting in 5nught/nught ^ j{AeughtY - (where /i is Hheavy)- 

These changes to the period accumulate over an entire cycle, giving 

Stiight ^ /i'e"'P. (2.31) 

By conservation of energy, the fractional change in semi-major axis (or period) of the 
heavy planet is reduced by a factor of mughtlf^heavyi so that 

^theavy - {might / mheavy) fJ'^ P- (2.32) 

2.4-2 Eccentricity dominated variations 

For timing variations that are caused by changes to the eccentricity of the transiting 
planet we look at the eccentricity dominated term in (j2.3(jp . Following the derivation 
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from the previous section (up to the point where the tisserand relation is apphed) we 
find that the eccentricity term gives a timing variation of 

St ~ /iperte-'P. (2.33) 

This means that for e > /i^/^ the perturbed eccentricity dominates the timing devi- 
ations while closer to resonance for j^/^jj?'/^ < e < ii^^"^ the perturbed mean motion 
dominates (this range is the same for both the light and heavy planets, except for fac- 
tors of order unity). For smaller values of e, the planets are trapped in mean- motion 
resonance, which was discussed in the previous section. Half way between resonances, 
e — so the timing deviation become 

5t ~ 0.7/ipert(a/(a2 - a^)fP. (2.34) 

A more precise derivation in the eccentricity-dominated case using perturbation the- 
ory is given in the next section. 

2.4-3 Perturbation theory derivation of eccentricity dominated variations 

In this section we consider more carefully the case of two planets whose orbits are 
nearly circular and whose timing variations are dominated by changes in eccentricity 
(see ^2.4|) . The timing variations can be computed from a Hamiltonian as described 
in j^sl. We keep terms which are first order in the eccentricity because mutual per- 
turbations between the planets induces an eccentricity of order m2/mQ. To first order 
in the eccentricities, the Taylor-expanded Hamiltonian is^ 

Hint= [2P-2«cos7/^-ei4cos(Ai-ti7i) 

+ 62^0 COS (A2 - ZU2) 

- ei J^JLi i^t " cos ((j - l)Ai - jA2 + uji) 

^We have corrected equation (26) in [s^ which should have a — acos')/' in the second line. 
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+ 62 E dj COS (jAi - (j - 1)A2 - rjj2) 

+ 62 E," 1 (4 - 4«(5,,i) COS (jAi - (j + 1)A2 + , (2.35) 

where 6ij is the Kronecker deha, cf = dab^y2 ^ "^3^^1/2' ~ '^f ~^ ^i/2('^)' ~ cli/(^2, 
by^i'^) is the Laplace coefficient, 



TT Jo vl^-2acosT+Q? 



h^],{a) = ^l de ,^ ^ , (2.36) 



where = Ai — A2, 



P = P(^/',a) = (l-2acos?/' + a^) (2.37) 



and 



dj^l{a)^a-^bf/,{a). (2.38) 

This equation includes no secular terms since these are higher order in the eccentric- 
ity. Note that since we have only included the ffist order terms in the eccentricity, 
the resonant arguments which appear have ratios j + l:j and j:j + 1 for the mean 
longitudes. 



39|, and we compute the perturbed 



The perturbed semi-major axis is given in 
eccentricity and longitude of periastron using hi = eisinroi and ki = eicostz?^. 
Keeping all the resonance terms that exist to ffist order in the eccentricities gives the 
equations of motion for hi,ki, 

h 1 "^2 r + . 

hi = — riia — CqCosAi 
2 mo 

00 

+ Yl { (4 - cos [(j - l)Ai - JA2] 
i=i 

+ {cj + adj^i) cos [(j + l)Ai - jA2] }] , 

ki = -Uia — [CgSmAi 

00 

5Z { (4 ~ sin [(j - l)Ai - JA2] 

j=i 

(cJ + 1) sin [(j + l)Ai - JA2]}] . 
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(2.39) 

To find the change in the transit timing we use the orbital elements to compute 
the variation in 6i. To first order in ei 

— ^ IH ^ + 2A;i cos [riiot + Aiol + 2hi sin [mot + Aiol. (2.40) 

nio nio 

Since we begin with zero eccentricity, we ignore perturbations to A in the sin and cos 
terms in this equation. As in equation ()2.23|) where 56*1 = 6i ~ riiQ, we integrate this 
equation to find 

Sti= 3^a-^^ (q(^) -asin^ 

^ mo (rai-n2) \^^^^'' 7r(l+a) ^ 

^ sin (Ai - Aio) 

I ^^"'^ ^^"■^'^'^J.i f sin [j'^/j] sin [nit+ji/)o] ^ 



^ mo 
-3a6 

j=l {j-l)ni-jn2 \j{ni-n2) ni J 

Eoo Cj +a5j,i / sin [jtp] _ sin [nit-jtpp] \ /2 ^j^N 

i=l (j+l)ni— jn2 "-2) ni y ' V • / 

where ni and n2 are taken at their initial values, Aio = \i(t = 0), A20 = A2(t = 0), 
k = 2y/a/{l + a), K{k) is the complete elliptic integral, Q{ip) is defined in the 
appendix of (sof, and we have dropped any terms which vary linearly with time. 

A similar calculation can be carried out for perturbations by a planet interior to 
the transiting planet, 



5t 



2 



mo (ni-?i2)' 

+ 

^ mo 



-«sin(^)) 



% sin[A2 - A20] 

'"2 

I Y^co lij / sm[jip] sin[n2f+jV'o] A 

Z^i=l {jni-{j+l)n2) \j{ni-n2) n2 J 

'^j 1 sin[jil>\ sin[n2<- 

{3ni-{j—l)n2) \j{ni-n2) n2 

A comparison of these equations to numerical calculations is shown in Figure The 
planets are on initially circular orbits with a semi-major axis ratio of 1.8, or period 
ratio of 2.4. The masses are equal and the planets start aligned along the line of sight 
at the first transit. 



_ spco d. / sin[jV»] _ s\a[n2t-j^l)o\ \ ^2 42) 
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Figure 2.3: Deviations from uniform times between transits for a transiting outer 
planet with semi-major axis ratio 1.8. Crosses represent the analytic result, equation 
()2.42|1 . while the diamonds are from a numerical integration of the equations of motion. 
The horizontal axis is the number of the transit, while the vertical axis shows the 
timing differences in units of Pi, the period of the inner perturbing planet, multiplied 
by mo/nipert- 
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In the case that the semi-major axis dominates the timing variations, one can take 
perturbed value of the eccentricity {h and k) and compute the change in semi-major 
axis due to these eccentricities. The result is a double-series over resonant terms, so 
I do not reproduce it here. 

So far we have discussed the timing variations for planets near, but not in a first 
order resonance. For larger period ratios, the eccentricity of the inner planet grows 
to Bin ~ Hout{Pin/Pouty, SO 6tin ~ IJ'outi.Pin / PoutY Pin- For an outer transiting planet 
the motion of the star dominates over the perturbation due to the inner planet for 
Pout > (27r)3/4p,„. 

Figure 12.41 shows a numerical calculation of the standard deviation of the transit 
timing variations. I used small masses to avoid chaotic behavior since resonant overlap 
occurs for j > yU"^/'^ Figure 1231 zooms in on the 2:1 resonance. As predicted, 

the amplitude scales as (equation I2.33|l . and then steepens to e~'^ (equations 12.321 
and 12.311) closer to resonance. Since the strength of the perturbation is independent 
of whether the perturbing planet is interior or exterior, the strength of the resonances 
are similar and the shape of the standard deviation of the transit timing variations 
is symmetric about Pj„ = Pout- The dashed curve in Figure 1231 shows the analytic 
approximation from equation ()2.10|) . which agrees well for Ppert < i'27r)^^/^Ptrans- The 
numerical results match the perturbation calculation, equations 1)2.411) and ()2.42|) . 
except for near resonance where the change in mean-motion dominates (we have not 
bothered to overplot the perturbation calculation since it is indistinguishable from 
the numerical results). 

There is a dip in a2 near Pout = 2.5Pj„ which occurs because the amplitude 
of the timing differences due to the orbit of the star about the barycenter (eqn. 
I2.10|) are opposite in sign and comparable in amplitude to the differences due to the 
perturbation of the outer planet by the inner planet (eqn. 12.42)) . 
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Figure 2.4: Transit timing standard deviation for two planets of mass rritrans = IQ-^mo 
and mpert = IQ-^mo on initially circular orbits in units of the period of the transiting 
planet. The solid (dotted) line is the numerical calculation for the inner (outer) planet 
averaged over 100 orbits of the outer planet with the planets initially aligned with 
the observer. The dashed line is equation ()2.10|) . The large dots are equations ()2.43|) 
on resonance and equation ()2.34|) halfway between resonances. 
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2.5 Planets in resonance 



The results of the previous section break down near each mean-motion resonance 
because hnear perturbation theory, which assumes that the perturbations are small, 
fails in the resonant regime. We must consider the effects of higher order terms and 
changes to the orbital elements of the perturbing planet in order to understand the 
effects of resonance. I present the case of low, initially zero, eccentricity where the 

r n 

standard analyses of this case [e.g. 49] to be incorrect. Here I provide a physically 
motivated, order of magnitude, derivation of the perturbations and the transit tim- 
ing variations for two planets in a first-order mean-motion resonance. A rigorous 
derivation is left for elsewhere, but I verify these findings with numerical simulations. 

Consider a first order, resonance where the lighter planet is a test particle. 

Qualitatively, the physics of low eccentricity resonance is as follows: on the nominal 
resonance, the two planets have successive conjunctions at exactly the same longitude 
in inertial space. The strong interactions that occur at conjunctions build up the 
eccentricity of the test particle and cause a change in semimajor axis and period. 
The change in period of the test particle causes the longitude of conjunction to drift. 
Once the longitude of conjunction shifts by about vr relative to the original direction, 
the eccentricity begins to decrease making a libration cycle. The libration of the 
semi-major axes causes the timing of the transits to change. Note that the libration 
of the longitude of conjunction distinguishes the resonant interaction with the non- 
resonant cases presented earlier. In non-resonant systems the longitude of conjunction 
constantly drifts in the same direction. 

The above, qualitative discussion leads directly to an estimate of the drifts in 
transit times. Within each libration cycle the longitude of conjunction shifts by about 
half an orbit, mostly due to the period change of the lighter planet. Since conjunctions 
occur only once every j orbits the largest transit time deviation of the lighter planet 
during the period of libration is P/j (in this order of magnitude derivation we ignore 
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factors of order unity, and take the limit of large j so that j ^ j + 1 and P2 — Pi). The 
observationally more interesting case is probably that in which the heavier planet is 
the transiting one. Then, conservation of energy for the orbiting planets implies that 
the change in periods is inversely proportional to the masses, therefore the timing 
variations are given by {rriught / fnheavyjP / j ■ We find an excellent fit to the data for 

Strr^a. ~ • (2.43) 

The calculations shown in Figure |231 verify this analytic scaling with j. 

Calculating the libration period is a little more complicated, but still straightfor- 
ward. Suppose the period of the test particle deviates from the nominal resonance 
by a small fraction e. Then, consecutive conjunctions drift in longitude by about 
27rj^e. The number of conjunctions, Nc, before a drift of order vr in the longitude of 
conjunctions accumulates is Nc ~ j~^e~^. We now estimate e indirectly. The test 
particle gains an eccentricity of order j^/i in each conjunction due to the radial force 
from the massive planet (this can be computed from the impulse approximation and 
the perturbation equation for eccentricity). The eccentricity given in N^. conjunctions 
is then of order Ae ~ yue^^. Using the Tisserand relation, the fractional change in 
semimajor axis associated with this change in eccentricity is j/i^e"^. Since this is also 
the fractional change in the period we have e ~ j^/^/i^/^ and a libration period of 

Pm ~ 0.5j-^e-ip ~ 0.5j-^/V"^^^^- (2-44) 

We numerically computed the amplitude and period of the transit timing varia- 
tions at the 2:1 resonance. Figure EUl shows a plot of the amplitude of the timing 
variations versus the mass ratio of the perturbing planet to the transiting planet. 
As predicted, the amplitude is of order the period of the transiting planet when the 
transiting planet is lighter, and varies as the mass ratio when the transiting planet 
is heavier. The libration period measured from the numerical simulations shows the 
predicted behavior, scaling precisely as ji'"^/^ for the more massive planet (with a co- 
efficient of ~ 0.7 for j = 1 and 0.5 for j > 1 in equation I2.44p . We have compared the 
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numerical values of the amplitude and period of libration on resonance as a function 
of j. Despite the fact that the above scalings were derived in the large-j limit, the 
agreement is better than 10 per cent for j > 2, and accurate to about 40 per cent for 

Figure 12.51 shows the more detailed behavior of the amplitude of the dispersion 
of the timing deviations near the 2:1 resonance. The amplitude is maximum slightly 
below resonance at the location of the cusp. This may be understood as follows: 
since the simulations are started with ei = 62 = 0, after conjunction the eccentricity 
grows and the outer planet moves outwards, while the inner planet moves inward. 
This causes the planets to move closer to resonance, causing a longer time between 
conjunctions, leading to a larger change in eccentricity and semi-major axis. The 
cusp is the location where the planets reach exact resonance at the turning point of 
libration, at which point 6t is maximum. To the right of the cusp, the libration causes 
the planets to overshoot the resonance, so the change in eccentricity and semi-major 
axis is somewhat smaller, and hence the amplitude is smaller. Figure ITHl shows that 
the width of the resonance scales as fi^^^ (the horizontal axis has been scaled with 
so that the curves overlap), so for larger mass planets the resonant variations have 
a wider range of influence than the non-resonant variations discussed in the previous 
section. The curves in Figure EiH demonstrate that on-resonance the amplitude scales 
as rnin{l, iipert/ l^trans) / i ^ while off-resonance the amplitude scales as /ipert- 

2.6 Non-zero eccentricities 

When either eccentricity is large enough, higher order resonances become important. 
In particular, the resonances that are l:m begin to dominate as the ratio of the semi- 
major axes becomes large; as the eccentricity of the outer planet approaches unity 
these resonances become as strong as first order resonances |5l| . Figure 12.61 shows 
the results of a numerical calculation where the transiting planet HD209458b, with 
a mass of approximately 0.67 Jupiter masses, is perturbed by a IM^ planet with 
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Figure 2.5: Amplitude near the 2:1 resonance versus the difference in period from 
exact resonance for two systems: one with mi = 10~^mo (top solid) and m2 = 10~^ 
(lower solid), and the other with mi = 10~^mo (top dotted) and m2 = 10"^ (lower 
dotted). The large dots are equation ()2.43|) . The dashed line is equation ()2.33|) . while 
the dash-dot line is equation ()2.32|) . 
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various eccentricities (we have taken HD209458b to have a circular orbit). Near the 
mean-motion resonances the signal is large enough that an earth-mass planet would 
be detectable with current technology. The amplitude increases everywhere with 
eccentricity. This graph can be applied to systems with other masses and periods 
as the timing variation scales as 5t oc PtransfT^pert (except for planets trapped in 
resonance). 

When both planets have non-zero eccentricity, the parameter space becomes quite 
large: the 4 phase space coordinates for each planet, assuming both are edge-on and 
2 mass-ratios give 10 free parameters. On resonance, the analysis remains similar 
to the circular case. The libration amplitude will still be of order ~ Pj^^ for the 
lighter planet and ~ Pj~'^{fi'iight/ fi'heavy) for the heavier planet. However, the period of 
libration will decrease significantly as the eccentricity increases since Puh — e~^/^yU~^''^. 

On the secular time-scale, the precession of the orbits will lead to a significant 
variation in the transit timing 48]. The period of precession, Pi^, may be driven by 
other planets, by general relativistic effects, or by a non-spherical stellar potential, but 
leads to a magnitude of transit timing deviation which just depends on the eccentricity 
for P,^ ^ P. 48] showed that the maximum deviation for e <^ 1 is given by 

6t=— (2.45) 

TT 

and the timing variations vary with a period that is equal to the period of precession. 
For arbitrary eccentricity, the maximum deviation is 

P 



St 

27r 



sin ""^y + sin ^ z + V2x — x"^ — x"^ , (2.46) 

where x = (1 — e^Y^^, y = (1 — x)/e, and z = {1 — x^)/e (this is derived from 
the Keplerian solution with a slowly varying w). This approaches P/2 as e ^ 1. 
Typically the eccentricity will vary on the secular time-scale, so these expressions 
only hold as long as the variation in e is much smaller than its mean value. 

Rather than systematically studying the entire parameter space, we now investi- 
gate several specific cases of known extrasolar planets to demonstrate that detection 
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Figure 2.6: Dispersion, atransi of timing variations of HD209458b due to perturbations 
induced by an exterior earth-mass planet with various eccentricities and periods. The 
color bar has (Jtrans for a planet of mass ttiq, and for a transiting period of 3.5 days. 
Near the period ratio of 4:3 the system becomes chaotic. The increase in the number 
of resonance peaks, particularly those beyond the 2:1 resonance, are from higher order 
terms in the expansion of the Hamiltonian which were truncated for the near-circular 
case. 
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Table 2.1: Timing variations for known multi-planet systems 



System 


P^n (d) 


Pout/ Pin 


0-1 


0-2 


55 Cnc e, b 


2.81 


5.21 


10.5 s 


2.68 s 


55 Cnc b, c 


14.7 


3.02 


1.61 h 


14.7 h 


Ups And b, c 


4.62 


52.3 


1.30 s 


1.61 min 


Gliese 876 


30.1 


2.027 


1.87 d 


14.6 h 


HD 74156 


51.6 


39.2 


4.98 min 


42.4 min 


HD 168443 


58.1 


29.9 


12.9 min 


2.62 h 


HD 37124 


152 


9.81 


3.43 d 


11.2 d 


HD 82943 


222 


2.00 


34.9 d 


30.7 d 


PSR 1257+12 b, c 


66.5 


1.48 


15.2 min 


22.3 min 


Earth/ Jupiter 


365 


11.9 


1.42 min 


24.1 s 



of this effect should be possible once a transiting multi-planet system is found. Most 
of these systems have non-zero eccentricities and several are in resonance, causing a 
significant signal. We summarize the amplitude of the signals of most known multi- 
planet systems, if they were seen edge-on, in Table 1 (in some cases other planets are 
present which would cause additional perturbations). 

The extrasolar planetary system Gliese 876 contains two Jupiter-mass planets on 
modestly eccentric orbits which are near the 2:1 mean-motion resonance. Pi = 30. Id 

r~i 

and P2 = 61. Od [42]. Due to the small size of the M4 host star, the inner planet has 
a 1.5 per cent probability of transiting for an observer at arbitrary inclination. The 
orbital motion involves both mean-motion resonance as well as a secular resonance 

alignment. The apsidal alignment is 
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50 



59I l37| . Figure 12.71 shows the 



in which the planets librate about their apsida 
in turn precessing at a rate of —41° per year 
predicted timing variations if this system were seen edge-on and if the planets are 
coplanar using the orbital elements from 
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The two most prominent periodicities in Figure EZI are those associated with the 
2:1 hbration, with a period of roughly 600 days [20 orbits of the inner planet, bd], and 
the long term precession of the apsidal angle with a period of about 3200 days (110 
orbits of the inner planet, corresponding to —41° yr^^). Evaluating equation ()2.45|1 
gives a peak amplitude of 1.4 days for the inner planet and 18 hours for the outer 
planet which both compare well with the numerical results given that the eccentricities 
are not constant. 

The extrasolar planetary system 55 Cancri contains a set of planets, b and c, near 
the 3:1 resonance having 15 and 45 day periods. There is some evidence for another 
planet, d, in an extremely long orbit, and recently a fourth low mass planet, e, was 
found with a 2.8 day period j46|. The planets e, b, and c have transit probabilities of 
12, 4, and 2 per cent, respectively, for an observer at arbitrary inclination. The orbit 
of planet b is approximately circular while planet c is somewhat eccentric [4^. Table 
1 gives the amplitude of the variations for the planets. We have ignored planet e; 
however, it is at a large enough semi-major axis to produce a ~ 22 second variation 
due to light-travel time as the barycenter of the inner binary orbits the barycenter of 
the triple system were the inner planets transiting. 

The double planet system Upsilon Andromedae has a semi-major axis ratio of 14 
which is not in a mean-motion resonance 16|, |42| . The inner planet has a short period 



of 4.6 days, and thus a significant probability of transiting of about 12 per cent, but 
has variations which are too small to currently be detected from the ground or space. 
The outer planet has much larger transit timing variations due to its smaller velocity, 
but a much smaller probability of transiting. 

The planetary system HD 37124 has two planets with a period ratio of ~ 10 and 
a period of the inner planet of 241 days ^69^. The outer planet is highly eccentric, 
62 = 0.69, and so its periapse passage produces a large and rapid change in the transit 
timing of the inner planet. If this system were transiting, the variations would be 
large enough to be detected from the ground. HD 82943 is in a 2:1 resonance giving 
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Figure 2.7: (a) Transit timing variations for Gliese 876 B and (b) Gliese 876 C. The 
vertical axis is in units of hours, while the horizontal axis is in units of the period of 
the transiting planet (given in parentheses for reference). 
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variations of order the periods of the planets. The pulsar planets are near a 3:2 
resonance, which would cause large transit timing variations were they seen to transit 
the pulsar progenitor star. Finally, alien civilizations observing transits of the Sun by 
Jupiter would have to have ~ 10 second precision to detect the effect of the Earth. 

2.7 Fourier Search for Semimajor Axis Ratio 

Since the TTV signal is typically periodic, it may be possible to use a discrete Fourier 
transform (FT) of the timing deviations to identify several of the orbital elements. An 
analysis conducted with a Fourier representation of the data may even be more appro- 
priate for finding the semimajor axis ratio of the system than analyzing the data in 
the time basis because the orbit of the periodic nature of the systems involved. Other 
orbital elements may also be determined with the Fourier representation, though such 
a development is left for elsewhere. 

Here, I present a method to determine the semimajor axis ratio of the planetary 
system using a Fourier transform of the timing deviations. The value of this tech- 
nique is that it could reduce the parameter space by one parameter entirely or, at 
least, provide several individual values about which a more complete search can be 
conducted. This technique may also give an independent confirmation of an estimate 
for the semimajor axis ratio parameter that is obtained by another means. 

2.7.1 Fourier Representation 

Consider figures ITHl and Figure ITHl shows the timing residuals of the transits of a 
planet where the perturber is near a mean-motion orbital resonance (the 2:1 resonance 
in this case). The second image in this picture shows the FT of those residuals. The 
symmetry in the Fourier transform is due to the fact that the timing information is 
encapsulated in both the amplitude and the phase of the Fourier components so that 
only the Fourier components up to the Nyquist frequency (half-way across the graph) 
contain unique information. For this case the Fourier component that corresponds to 
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the 2:1 resonant forcing term has the largest amphtude. By comparison, figure EHl is 
for a system that is not near an orbital resonance and has several large peaks in the 
amplitudes of the Fourier components. 

In order to determine whether or not the relative spacing between the peaks in 
the FT was interesting, I generated a set of planetary systems where both planets 
had zero eccentricity, equal mass, equal time of pericenter passage, and with various, 
equally spaced semimajor axis ratios with the perturber being exterior to the transit- 
ing planet. For each system I tabulated the variations in transit time and then took 
the FT of the timing residuals. I then calculated the time intervals that correspond 
to the largest peaks (up to five) of the FT. 

From these results I generated a scatter plot of the values of the five time intervals 
as a function of the semimajor axis ratio and obtained figure EH Each of the visible 
peaks corresponds to a particular mean-motion resonance — the shape of which is 
given later. One of the remaining structures is a diagonal line that corresponds to the 
period of the perturbing planet in units of the period of the transiting planet. This 
diagonal line is particularly interesting because it is invertable and can, therefore, be 
used to uniquely identify the period of the perturbing planet. 

The overall structure of this graph is virtually independent of the orbital elements 
of the system. For example, another set of systems, with identical values of the semi- 
major axis ratio but with the remaining orbital elements being randomly generated 
results in Figure 12.111 The primary difference between Figures I2.1UI and 12.111 is that 
the relative heights of the various peaks in the FT; the peaks themselves remain in 
the same locations. 

The shape of the peaks that are visible in these two plots can be derived from 
the corresponding mean motion resonance. For the r.j resonance, where i and j are 
integers and where iPtrans = j-Ppcrt corresponds to exact resonance, the period that 
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Figure 2.8: Timing residuals (upper panel) with Fourier transform of those residuals 
(lower panel) for a planetary system that is near a mean-motion resonance. 
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Figure 2.9: Timing residuals (upper panel) with Fourier transform of those residuals 
(lower panel) for a planetary system that is not near a mean-motion resonance. 
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Figure 2.10: This picture shows the periods that correspond to the five largest peaks 
in the Fourier transform of the timing residuals for a transiting planet where both 
planets are on initially circular orbits as a function of the ratio of the semimajor axes 
of the planets. The period that corresponds to the highest peak is in red, the second 
highest is in yellow, then green, blue, and purple. 



39 



O 1 

0.5 



'• 

• • 








• 




• • • • 






• • • 




• • 


• • • 


• 


• 


■ • 




• 


• 


• 




• 

• • • 


• 


• 


• • • • • • 

• • • ••• • 


• 


• 


• • • • • • • 


• 

• < 

1 


• 




















• • • 
.* ••• • 

• ••• . • • 


li 


:i r*'*' ::: 















1.5 1.75 2 2.25 2.5 2.75 3 
Axis Ratio 



Figure 2.11: Similar to Figure l2.10| this shows the periods that correspond to the 
five largest peaks in the Fourier transform of the timing residuals as a function of 
the semimajor axis ratio of the two planets. Here the orbital elements of the systems 
at each axis ratio are randomly selected. Red corresponds to the highest peak, then 
yellow, green, blue, and purple. 
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Figure 2.12: This figure has the plot from Figure 12.101 with several curves of equation 
12.471 superposed. Also shown are the two diagonal lines that correspond to the period 
of the perturbing planet (upper diagonal line) and half the period of the perturbing 
planet (lower diagonal line). 



describes the shape of the structures in Figures 12.101 and 12.111 is given by 

-Ppcak _ -Ppcrt ^2 ^-^^ 

-Ptrans ^-Ptrans jPpevt 

This shows that the period of the Fourier component that appears on a given peak 
approaches infinity as the system approaches the resonance and falls away as the 
distance from that resonance increases. Figure 12.121 superposes several curves of 
equation (|2.47|) . that correspond to different resonances, onto Figure 12. 101 

The value of this development is that the Fourier information can be used as a 
fingerprint to identify the period of the perturbing planet directly from the transform 
of the timing residuals. Starting with the FT of a given set of timing residuals we 
plot the periods that are associated with the peaks in the amplitude of the FT. We 
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Figure 2.13: This figure has the plot from Figure l?. lOl with the periods that correspond 
to the highest peaks in the Fourier transform — shown in black — of the resonant system 
whose transit times are shown in Figure |TH1 (on the left) and the nonresonant system 
whose transit times are shown in Figure 12 .91 



then assume that one of the peaks is precisely the period of the perturbing planet 
and compare the remaining peaks with the set of periods that result from equation 
(|2.47|) . For example, in Figure E. 131 1 show the periods that correspond to the peaks 
in the FT of the systems shown in Figures 12.81 and 12.91 superposed upon Figure I2.1UI 

2.7.2 Implementation 

Several tasks remain to be done before this technique to determine the semimajor 
axis ratio of the planetary system can be implemented in a robust and automated 
manner. The largest impedement is that I have not found an appropriate statistical 
test to estimate the proper fit of the Fourier peaks to equation (I2.47j) . I initially tried 
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several systems "by eye" and found that I could correctly identify approximately 90% 
of the systems — though no noise had been added to the transit timings. Most of the 
remaining systems had several possibilities that looked as though they fit equally well. 

The difficulty lies in the fact that there are infinitely many resonances against 
which the data could be tested; but they are not equally important. One could 
probably find a perfect fit to the Fourier peaks if he were allowed to include the effects 
of resonances such as the 6:3 resonance — which corresponds to the 2:1 resonance but 
still provides an independent curve from equation ()2.47j) (this can be seen in Figure 
I2.1()|l — or the 200:101 resonance which is also near the 2:1 resonance but which is 
intuitively less important. There is no clear method to determine which resonances 
should be weighted more heavily than others. One potentially fruitful method would 
be to weight each higher order resonance (where the order is defined by \i — j\) by 
the maximum, average, or product of the eccentricities of the planets. Thus, for low- 
eccentricity systems the importance of the high-order resonances would be severely 
limited; while for highly eccentric systems they would be more important. 

Another challenge is that there are not always multiple peaks that can be identified 
in the FT — see the Figures [2.81 and Any noise that is added to the data will tend 
to wash out one or more of the peaks regardless of its importance to determining 
the semimajor axis ratio parameter. For example, in the resonant system of Figure 
12.81 the second largest peak corresponds to the (unique) period of the perturbing 
planet. This peak is significantly smaller than the peak that corresponds to the 2:1 
resonance; but without it there is no way to identify that the resonance is indeed 
2:1 — any resonance would be a candidate. I have yet to determine an appropriate 
way to handle this issue and the one discussed in the preceeding paragraph. Until 
such time as these problems have been solved the value of this Fourier technique is 
uncertain. It does, however, illustrate the potential that the technique posseses. 



43 



Chapter 3 



APPLICATIONS AND CAVEATS 



3.1 Detection of terrestrial planets 

The possibility of detecting terrestrial planets using the transit timing technique 
clearly depends strongly on (1) the period of the transiting planet; (2) the proximity 
to resonance of the two planets; (3) the eccentricities of the planets. The detectabil- 
ity of such planets also depends on the measurement error, the intrinsic noise due to 
stellar variability, gaps in the observations, etc. One requirement for the case of an 
external perturbing planet is that observations should be made over a time longer 
than the period of the timing variations, which can be longer than the period of the 
perturbing planet. Ignoring these complications, a rough estimate of detectability 
can be obtained from comparing the standard deviation of the transit timing with 
the measurement error. 

It is worthwhile to provide a numerical example for the case of a hot Jupiter with 
a 3 day period that is perturbed by a lighter, exterior planet on a circular orbit in 
exact 2:1 resonance. The timing deviation amplitude is of order the period (3 days) 
times the mass ratio (300) or about 3 minutes fequation 12. 43^ : 



These variations accumulate over a time-scale of order the period (3 days) times the 
planet to star mass ratio to the power of 2/3, which for a transiting planet of order a 
Jupiter mass is about 5 months fequation 12. 44^ : 




(3.1) 




(3.2) 
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Such a large signal should easily be detectable from the ground. With relative pho- 
tometric precision of 10~^ from space or from future ground-based experiments, less 
massive objects or objects further away from resonance could be detected. The ob- 
servations could be scheduled in advance and require a modest amount of observing 
time with the possible payoff of being able to detect a terrestrial-sized planet. 

3.1.1 Comparison to other terrestrial planet search techniques 

To attempt a comparison with other planet detection techniques, we have estimated 
the mass of a planet which may be detected at an amplitude of 10 times the noise for 
a given technique. We compare three techniques for measuring the mass of planets: 
(1) radial velocity variations of the star; (2) astrometric measurements; (3) transit 
timing variations (TTV). We assume that radial velocity measurements have a limit 
of 0.5 m/s RMS uncertainty, which is about the highest precision that has been 
achieved from the ground, and may be at the limit imposed by stellar variability 
We assume that astrometric measurements have a precision of 1 /larcsecond 
which is the precision which is projected to be achieved by the Space Inter ferometry 
Mission [g^ . Finally, we assume that the transit timing can be measured to a 
seconds, which is the highest precision of transit timing measurements 



precision of IC 
of HD209458 



3- 

We concentrate on HD209458 since it is the best studied transiting planet. This 
system is at a distance of 46 pc and has a period of 3.5 days. Figure 13.11 shows a 
comparison of the sensitivity of these three techniques for a signal to noise ratio of 
10. The solid curve is computed for rritrans = 6.7 x IO^^Mq and mp^rt = 1O^''M0 
and both planets on circular orbits. When not trapped in resonance, the amplitude 
of the timing variations scales as nipert/rno, so we scale the results to the mass of the 
perturber to compute where the timing variations are 100 seconds - this determines 
the sensitivity. The TTV technique is more sensitive than the astrometric technique 
at semi-major axis ratios smaller than about 2. Note that in Figure IXTl the TTV and 
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Figure 3.1: Greatest lower bound on the detectable mass of a secondary planet for 
various planet detection techniques for the HD209458 system. The vertical axis is 
the perturbing planet's mass in units of Mq. The horizontal axis is the period ratio 
of the planets. The solid line is for the transit timing technique, the dashed line is 
astrometric, and the dotted line is the radial velocity technique. 



astrometric techniques have the same slope at small Ppert / Ptrans ■ This is because the 
transit timing technique is measuring the reflex motion of the host star due to the 
inner planet, which is also being measured by astrometry and that the solid curve is an 
upper limit to the minimum mass detectable in HD209458 since a non-zero eccentricity 
will lead to larger timing variations (Figure 12. 6|) and thus a smaller detectable mass. 

We see from l3.ll that. off resonance, radial velocity measurements are the technique 
of choice for this system, while on resonance the TTV is sensitive to much smaller 
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planet masses. The improvement in sensitivity with the TTV technique over the 
radial velocity technique for a typical j:j + 1 resonance is estimated as follows. The 
radial velocity that a planet induces in its host star scales as f ~ {a/ P){m„/mQ) 
where a is the semimajor axis, P is the period, and m^^ and mo are the planet 
and star masses respectively. By replacing the velocity of the star with the radial 
velocity measurement precision {v ^ a^) and solving for the mass of the planet, we 
see that the radial velocity technique is able to detect planets with masses of order 
m-rv ~ c^^,moP/a. 

Similarly, the typical timing deviations for a resonant system is given by equation 
12.431 with is approximately 5t ~ (P/10)(mttv/^p) where P is the period, nip is the 
transiting planet mass, and mttv is the perturbing planet mass which I assume is much 
smaller than mp. The factor of 10 comes from the denominator of equation 12.431 By 
replacing the transit timing deviation by the timing measurement uncertainty 5t — > cr^ 
and solving for the mass of the perturber we obtain the detection limit for the TTV 
technique mttv ~ lOatnip/ P. By dividing the mass limit for TTV with the mass limit 
for RV we find 

mttv lOcTta ( mp\ ,f_^^^rx'^p 



(few)io:::^ (3.3) 
m^v f^-T J rriQ 

for typical values of the different observational and transiting hot j up iter parameters. 
The large improvement (an order-of-magnitude or more) stems from the fact that the 
resonant TTV signal scales as the mass ratio of the planets because of conservation 
of the energy of the planetary orbits while the radial velocity scales with the planet 
to star mass ratio due to conservation of the barycenter motion. 

3.2 Breaking the mass-radius degeneracy 

In the case that two planets are discovered to transit their host star, measurement 
of the transit timing variations can break the degeneracy between mass and radius 
needed to derive the physical parameters of the planetary system. This has been 
discussed by jfiSj who use a theoretical stellar mass-radius relation to break this 
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degeneracy. We provide a simplified treatment here to illustrate the nature of the 
degeneracy and how it can be broken with observations of transit timing variations. 

Consider a planetary system with two transiting planets on circular orbits which 
are coplanar, exactly edge-on, and have measured radial velocity amplitudes. We'll 
assume that the star is uniform in surface brightness and that mi,m2 ^ mo. We'll 
also assume that the unperturbed periods Pi,P2 can be measured from the dura- 
tion between transits. Then there are eight physical parameters of interest which 
describe the system: mo,mi,m2, Rq, Ri, R2,cii, and 02 where Ri are the radii of the 
star and planets. Without measuring the transit timing variations, there are a total 
of ten parameters which can be measured: Ki,K2, tTi,tT2, tgi, tg2, AFi, AF2, Pi, and 
P2, where tTj labels the duration of transit from mid- ingress to mid-egress, tgj labels 
the duration of ingress or egress for planet j, Kj are the velocity amplitudes of the 
two planets, and AFj are the relative depths of the transits in units of the uneclipsed 
brightness of the star (for planet j). Although there are more constraining parame- 
ters than model parameters, there is a degeneracy since some of the observables are 
redundant. All of the system parameters can be expressed in terms of observables 
and the ratio of the mass to radius of the star, tjiq/Rq, 
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(3.4) 

where j = 1,2 labels each planet. From this information alone one can constrain the 
density of the star 6^ . For the simplified case discussed here. 
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for either planet [this differs shgthly from the expression in I63II, since we define the 
transit duration from mid in/egress]. If, in addition, one can measure the amphtude 
of the transit timing variations of the outer planet, (72, then this determines the mass 
ratio. For the case that the star's motion dominates the transit timing, 

nil 

For other cases, the transit timing amplitude can be computed numerically. Then, 
from the above expression for mj/mo one can find the ratio of the mass to the radius 
of the star 

Ro Stt-^G trial ' ^ ' ' 

Combined with the measurement of the density, this gives the absolute mass and 
radius of the star. This procedure requires no assumptions about the mass-radius 
relation for the host star, and in principle could be used to measure this relation. 
If one can also measure transit timing variations for the inner planet, then an extra 
constraint can be obtained 

a, = ^fia), (3.8) 
mo 

where /(a) is a function derived from averaging equation l2.41l (Note that the phase of 
the orbits is needed for this equation, which can be found from the velocity amplitude 
curve). This provides an extra constraint on the system, and thus will be a check 
that this procedure is robust. 

Clearly we have made some drastically simplifying assumptions which are not true 
for any physical transit. The inclination of the orbits must be solved for, which can 
be done from the ratio of the durations of the ingress and transit and the change 
in flux, as demonstrated by 6^. In addition, limb- darkening must be included, and 
can be solved for with high signal-to- noise data as demonstrated by Finally, 
the orbits are not generally circular, so the parameters ej,u7j,Qj, aj, which can be 
derived from the velocity amplitude measurements, should be accounted for. The 
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general solution is rather complicated and would best be accomplished numerically, 
but the degeneracy has a similar nature to the circular case and can in principle be 
broken by the transit timing variations. 

3.3 Effects we have ignored 

I now discuss several physical effects that I have ignored but which ought to be kept 
in mind by observers measuring transit timing variations. 

3.3.1 Light travel time 

j^l carried out a search for perturbing planets in the eclipsing binary stellar system, 
CM Draconis, using the changes in the times of the eclipse due to the light travel time 
to measure a tentative signal consistent with a Jupiter-mass planet at ~ 1 AU (their 
technique would in principle be sensitive to a planet on an eccentric orbit as well, c.f 
equation 12. 25|1 . The "R0mer Effect" due to the change in light travel time caused by 
the reflex motion of the inner binary is much smaller in planetary systems than in 
binary stars since their masses and semi-major axes are small, having an amplitude 

t^^^^ = P ^ O.Ssec — ( — ^ I 3.9 

cm, \MjJ VlAU/ ^ ' 

where Mj is the mass of Jupiter and 02 is the semi-major axis of the perturbing 
planet. This effect is present in the absence of deviations from a Keplerian orbit 
because the inner binary orbits about the center of mass. 

There can also be changes in the timing of the transit as the distance of the 
transiting planet from the star varies. In this case, the time of transit is delayed by 
the light travel time between the different locations where the planet intercepts the 
beam of light from the star. The amplitude of these variations is smaller than the a we 
have calculated by a factor of ~ Vtrans/c, where Vtrans is the velocity of the transiting 
planet. So, only very precise measurements will require taking into account light 
travel time effects, which should be borne in mind in future experiments (of course 
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the light-travel time due to the motion of the observer in our solar system must be 
taken into account with current experiments). 

3.3.2 Inclination 

We have assumed that the planets are strictly coplanar and exactly edge-on. The 
first assumption is based on the fact that the solar system is nearly coplanar and 
the theoretical prejudice that planets forming out of disks should be nearly coplanar. 

n 

Small non-coplanar effects will change our results slightly 48], while large inclination 
effects would require a reworking of the theory. Since some extrasolar planetary 
systems have been found with rather large eccentricities, it is entirely possible that 
non-coplanar systems will be found as well, a possibility we leave for future studies. 

The assumption that the systems are edge-on is based on the fact that a transit can 
occur only for systems that are nearly edge-on. For small inclinations our formulae 
will only change slightly, but may result in interesting effects such as a change in the 
duration of a transit, or even the disappearance of transits due to the motion of the 
star about the barycenter of the system. On a much longer time-scale (centuries), the 
precession of an eccentric orbit might cause the disappearance of transits since the 
projected shape of the orbit on the sky can change. This possibility was mentioned 
by y for GJ 876. 

3.3.3 Other sources of timing "noise" 

Aside from the long term effects that have been ignored there are several sources 
of timing noise that must be included in the analysis of observations of transiting 
systems. These sources of noise could come from the planet or the host star. If the 
planet has a moon or is a binary planet then there is some wobble in its position 
causing a change in both the timing and duration of a transit Il3l |. A moon or 



ring system may transit before the planet cau sing a shallower transit to appear earlier 



or later than it would without the moon 



13L l62l . |5[ . A large scale asymmetry of the 
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planet's shape with respect to its center of mass might cause a shghtly earher or later 
start to the ingress or end of egress. 

Stellar variability could also make a significant contribution to the noise. Varia- 
tions in the brightness of the star might affect the accuracy of the measurement of 
the start of ingress and the end of egress, which are the times that are critical to 
timing of a transit. Stellar oscillations can cause variations in the surface of the Sun 
of ~ 100 km in regions of size 10'^"'^ km, which corresponds to a one second variation 
for a planet moving at 100 km/s. 

3.3.4 Coverage Gaps 

With radial velocity measurements and prior transit lightcurves one can predict the 
epoch of future transits and identify appropriate times for photometric monitoring of 
the system of interest. Observational limitations (e.g. bad weather, equipment failure, 
scheduling requirements) will lead to transits being missed, which in turn will cause 
inaccuracies in a. Since the signal is periodic, St{t) may be straightforward to extract 
with a few missed transits; however, if the outer planet is highly eccentric, then most 
of the change in transit timing may occur for a few transits [e.g. HD37124; a similar 
selection effect occurs in radial velocity searches as discussed by|2o|. In principle this 
effect will average out over long observational intervals; however as in this context 
"long" may mean several decades or more, it will be important to evaluate the effect 
of coverage gaps on detections over a time-scale of months-years. We will return to 
this in detail in future work. 

We note in passing that the advent of the new astrometric all-sky surveys such 
as Gaia 54] will provide photometry for ~ 10"^ stars > 15 mag, with fewer cover- 
age gaps than ground-based observations; we thus expect the detection method by 
transits alone f section l3.1|) to really come into its own over the next two decades. As- 
suming ~ 0.4 detections (three transits) per 10^ stars [c.f. _14|, we may expect transit 
lightcurves of perhaps ~ 1000 exoplanetary systems over the mission lifetime, greatly 
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aiding the determination of a for these systems. The NASA Kepler mission will also 
provide uniform monitoring of about 10^"^ transiting gas giants Jo^, and if flown, the 
Microlensing Planet Finder (formerly GEST) will discover 10^~^ transiting gas giants 
with uniform coverage 

3.4 Conclusions 

For an exoplanetary system where one or more planets transit the host star the timing 
and duration of the transits can be used to derive several physical characteristics of 
the system. This technique breaks the degeneracy between the mass and radius of the 
objects in the system. The inclination, absolute mass, and absolute radii of the star 
and planets can be found; in principle this could be used to measure the mass-radius 
relation for stars that are not in eclipsing binaries. 

In addition, TTV can be used to infer the existence of previously undetected 
planets. We have found that for variations which occur over several orbital periods 
the strongest signals occur when the perturbing planet is either in a mean-motion 
resonance with the transiting planet or if the transiting planet has a long period 
(which, unfortunately, makes a transit less probable). The resonant case is more 
interesting since the probability of a planet transiting decreases significantly as the 
semi-major axis becomes large. Using the TTV scheme it is possible to detect earth- 
mass planets using current observational technologies for both ground based and space 
based observatories. Observations for several transits of HD209458 could be gathered 
and studied over a relatively short time due to its small period. Once the existence 
of a second planet is established one can predict the times at which it would likely 
transit the host star. Follow-up observations with HST or high-precision ground- 
based telescopes at those times would increase the likelihood of detecting a transit of 
the second planet. 

If the second planet is terrestrial in nature, this transit timing method may be 
the only way currently to determine the mass of such planets in other star systems. 
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Astrometry is another possible technique but it may take a decade of technological 
development before the necessary sensitivity is achieved. In addition, complementary 
techniques are necessary to probe different parts of parameter space and to provide 
extra confidence that the detected planets are real, given the likely low signal-to-noise 
For the near future the TTV technique may be the most promising method of 
detecting earth-mass planets around main sequence stars besides the Sun. 

We exhort observers to (1) discover longer period transiting planets 6^ since 
the signal increases with transiting planet's period; (2) increase the signal-to-noise of 
ground based differential photometry j^^] for more precise measurement of the transit 
times; and (3) examine their transit data for the presence of perturbing planets jisl . 

The treatment of this problem has ignored many effects which we plan to take 
into account in future works in which we will simulate realistic lightcurves including 
noise and to fit the simulated data to derive the parameters of the perturbing planet, 
exploring degeneracies in the period ratio. We will also derive the probability of 
detecting such systems taking into account various assumptions about the formation, 
evolution, and stability of extrasolar planets. At this stage we do not address the 
effects of general relativity since the timescale over which these effects would be 
manifest is significantly longer than the timescale over which existing and planned 
observations are to collect data. 
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Chapter 4 
THE TRES-1 SYSTEM 

In this chapter and the next I present two analyses of existing transit data. The 
first data set is for the TrES-1 planetary system (GSC 02652-01324) and the second 
is for the HD209458 planetary system. The transit data for the TrES-1 system were 
obtained using ground-based telescopes and were reported by David Charbonneau 
in I19I] (hereafter C05). The data for the HD209458 system were obtained with the 
Hubble Space Telescope. Many of the results of this chapter are reported in {g^. 

The timing data reported in C05 were derived from the 11 transits reported by J] 
with an additional transit that was observed at the lAC 80cm telescope after J] went 
to press. One transit was excluded from their analysis because it constituted a 6-a 
departure from a constant period and because of anomalous features in the ingress 
and egress. That point, if it is valid, is the most interesting point for our purposes 
because the TTV signal is defined by such deviations. Consequently, I analyse two 
different sets of data from C05; the "12-point" set which includes this point, and the 
"11-point" set which does not. This study was the first analysis of TTV as presented 
by 1| and j^]. And, it may be the first search for planets around main-sequence 
stars that can probe masses smaller than the mass of the Earth. 

4.1 Search for Secondary Planets 

I conducted a variety of searches for the best-fitting perturbing planet in the TrES-1 
system. These searches included different combinations of orbital elements for TrES- 
Ib. I found that any reduction in the overall obtained by including the parameters 
e and w for TrES-lb was offset by the loss of a degree of freedom. Therefore, I report 
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results from the search where the eccentricity of TrES-lb was fixed at zero. 

For this analysis, I stepped through the semi-major axis ratio of the putative 
secondary planet and TrES-lb. At each point I minimized over six parameters: the 
eccentricity, longitude of pericenter, time of pericenter passage, and mass of the sec- 
ondary planet and the period and the initial longitude of TrES-lb. The inclination 
and ascending node of the perturbing planet were identical to the values for TrES-lb. 
I analysed the data for both interior and exterior perturbers. 

This analysis did not produce any promising solution; though I present one inter- 
esting case for an interior perturber found from the 12-point analysis. This solution 
is not near a low-order mean-motion resonance. Indeed, the was generally much 
higher near the j:j + 1 and the j:l resonances than in the regions between them. 
Figure compares the timing residuals for this solution with the data. The reduced 

for this system is 2.8 on — 6 degrees of freedom (where iV = 12 is the number 
of data) compared with 6.3 for no perturber (A^ — 2). 

This solution, while it would have been detected from RV measurements, is in- 
teresting because the average size of the timing deviations is larger than that of the 
data — making the variations easier to detect. However, I suspect that while this solu- 
tion is numerically valid, it is merely an artifact of the gap between the two primary 
epochs of observation. For this solution, and others like it, the simulated timing 
residuals consist of small, short-term variations superposed on a large-amplitude, 
long-period variation with a period that is a multiple of the difference between the 
two epochs. Several candiate systems for both the 11-point and the 12-point analyses 
showed similar behavior. Such results indicate one drawback of gaps in the coverage 
of transiting systems. However, observing each transit of every transiting system may 
be neither feasable nor optimal for identifying perturbing planets with this technique 
given the limited resources that are appropriate for transit observations. In the mean- 
time, additional data for this system, taken at a time that is not commensurate with 
the existing gap in the observations would remove false solutions of this type in future 
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Figure 4.1: A comparison of the TrES-lb transit timing residuals to those from the 
12-point analysis for an interior perturber. The diamonds are the timing residual at 
each transit, the filled diamonds denote transits where there is data, the squares with 
error bars are the data. The 12-point solution has a reduced of 2.8 on 6 degrees of 
freedom with ephemeris parameters m — O.lllMj, P — 1.76d, e = 0.153, uj — 154.8°, 
T = 2452847.4392d. 



studies. 

Prom the above analysis, I conclude that there is not sufficient information in 
the data to uniquely and satisfactorally determine the characteristics of a secondary 
planet in the TrES-1 system. This is in part because the number of model parameters 
is not much larger than the number of data and because the typical timing precision, 
~ 100s, is not a sufficiently small fraction of the orbital period of the transiting 
planet (about 4 x 10~^) to distinguish between different solutions. In addition, the 
gap in coverage appears to affect the minimization dramatically. I believe that the 
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coverage gap is primarily responsible for the observed fact that two nearby sets of 
orbital elements will have very different values of — a slight change in the long- 
term variation will cause the simulated transit times to miss the second epoch of 
observations. Additional timing data, with precision comparable to the most precise 
of the given data, ~ 30s, and at an epoch that is not commensurate with the existing 
gap in coverage will allow for a more complete investigation of the system. 

4.2 Constraints on Secondary Planets 

That the results of the planet search were inconclusive is not particularly surprising 
since the transit timings show little variation from a constant period and the point 
that deviates the most is suspect. If a satisfactory solution cannot be determined from 
these data, we can still use them to place constraints on the orbital elements of a two- 
planet system. Of particular interest is a constraint on the mass that the secondary 
planet could have as a function of various orbital elements; similar to the constraint 
that resulted from analyzing the OGLE-1998-BUL-14 microlensing event (sj. As the 
mass of the hypothetical perturbing planet increases, the of the model should grow 
significantly regardless of the values of the remaining orbital elements. Therefore, I 
systematically studied a grid of values on the semi-major axis/eccentricity plane of 
the perturbing planet — assuming that the orbit of TrES-lb is initially circular. At 
each point we identified the mass that a perturber needs in order to produce a large 
deviation from the data. 

At each location in the a/e plane we set the mass of the companion to be very 
small (<^ Mq^) and, for a random value for the longitude of pericenter and the time 
of pericenter passage, we calculated the of the timing residuals. We increased the 
mass of the secondary planet until the grew by some fiducial amount. At that 
point we minimized the over the longitude of pericenter and the time of pericenter 
passage of the secondary planet and the initial longitude of TrES-lb. Following the 
minimization the typically fell below the fiducial amount and the mass of the 
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perturber was again increased. However, if after the minimization the remained 
above the threshold, the orbital elements of the system were recorded, the threshold 
increased, and the process repeated until either a maximum mass of the perturbing 
planet or a maximum was achieved. This procedure gives the minimum as a 
function of perturbing planet mass for each point in the a/e plane. 

The mass that corresponds to a 3-0" increase in the constitutes our estimate 
for the maximum allowed mass of a secondary planet. I use the approach outlined 
by where I locate the mass that causes the difference between the of the null 
hypothesis and the x^ obtained with a secondary planet to equal Ax^ = 9. This 
gives the S-cx limit because only one parameter is allowed to vary — the remaining 
parameters are either fixed or marginalized at each point. This maximum mass is 
shown as a function of a and e for the 11-point set in Figures f4. 2114. 31 for an interior 
and an exterior perturber respectively. The contours correspond to 100, 10, and 1 
Mq,. The dark portion in the upper- left corners are where the orbits overlap and I 
assign a mass of 1Q~^Mq to those locations. We see from these figures many regions 
where the mass of a companion must be comparable to or less than the mass of the 
Earth regardless of its orbital eccentricity. The most stringent constraints are near the 
j:j + 1 mean- motion resonances. Of particular interest are the very tight constraints 
for low-eccentricity perturbers (which constraints continue to an eccentricity of zero) 
because tidal circularization drives the eccentricity down. The 12-point analysis gave 
similar results which shown in Figures f4.4l and WJ\ 

4.3 Comparison with Radial Velocity 

Much of the region where a secondary planet is not tightly constrained by this analysis 
could be limited more efficiently by radial velocity measurements. Figure 14.61 shows 
the constraint achieved from this TTV analysis and the constraint from a hypothetical 
RV analysis, with the same number of data, as a function of the period ratio of the 
two planets. I assume that the RV measurements have a precision of both 5 and 1 
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Figure 4.2: The upper limit on the mass of an interior perturbing planet in the TrES-1 
system from the 11-point analysis. The regions bounded by the contours correspond 
to a companion mass that is less than 100 (dotted), 10 (dashed), and 1 (solid) earth 
masses. 
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Figure 4.3: The upper limit on the mass of an exterior perturbing planet in the TrES-1 
system from the 11-point analysis. The regions bounded by the contours correspond 
to a companion mass that is less than 100 (dotted), 10 (dashed), and 1 (solid) earth 
masses. 
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Figure 4.4: The upper limit on the mass of an interior perturbing planet in the TrES-1 
system from the 12-point analysis. The regions bounded by the contours correspond 
to a companion mass that is less than 100, 10, and 1 earth masses. 
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Figure 4.5: The upper limit on the mass of an exterior perturbing planet in the TrES-1 
system from the 12-point analysis. The regions bounded by the contours correspond 
to a companion mass that is less than 100, 10, and 1 earth masses. 
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m and that the orbit of the perturber can be treated as circular. Figures 14.71 
and 14. 81 are the same comparison focussed on the region surrounding the 2:1 and 3:2 
mean-motion resonances respectively. 

I derived the RV limit in these figures by assuming that the velocity residuals 
that remain after removing the effects of TrES-lb randomly sample the phase of the 
putative second planet and that the noise is uncorrelated with the RV measurements. 
If no secondary planet exists, then the residuals should surround zero with a variance 
equal to cr^ where a is the precision of the measurements. The expected is then 
equal to = N — Qq where Qo is the number of model parameters. On the other 
hand, if a secondary planet exists, the expected is 

2/3 

2 V ^ 

where /i is the planet to star mass ratio, Mq is the mass of the star, G is Newton's 
constant, and Q is the number of model parameters. By subtracting Xq from and 
solving for the mass ratio we obtain the mass of a secondary planet that is detectable 
with a given Ax^ threshold as a function of the period of the planet and the precision 



2 {N-Q) 

^2 



/i2 f2TlGMo \ ^ ^2 



(4.1) 



of the radial velocity measurements 



/i = O" 



2(Ax' + Q-gn)r'v p 



1/2 / p ^ 1/3 

(4.2) 



2nGMn 



N-Q 

We see from Figure 14.61 that for low-order, mean-motion resonances the TTV 
analysis is able to place constraints on the mass that are nearly an order-of-magnitude 
smaller than the RV technique with 5 m s~^ precision and that there are regions where 
it is more sensitive than RV measurements with 1 m s~^ precision. Meanwhile, 
in nonresonant regimes the RV method remains superior for a large portion of the 
parameter space. Additional transit timing data, particularly with precision that is 
comparable to the most precise of the given data, would lower the entire limit from 
TTV, provided that no secondary planet exists. Such data would render the TTV 
approach superior over a larger range of periods. 
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Figure 4.6: A comparison of the limits expressed in this paper for the 11-point analysis 
(solid curve) and the expected limits that are achieved with 11 randomly selected 
radial velocity measurements that have a precision of 5 m s~^ (dotted curve) and 1 
m s~^ (dashed curve) for a secondary planet with an eccentricity of 0.05. The radial 
velocity curves are given by equation ()4.2|) using Ax^ = 9, N = 11, Qo = 0, and 
Q = 2. Near the 4:3, 3:4, and the 1:2 mean-motion resonances the TTV analysis of 
these data can place tighter constraints (about a factor of 3) on the mass of putative 
secondary planets than can the RV technique with 1 m s~^ precision. 



65 



Interior Perturber 



- 3 




0.49 0.50 0.51 0.52 
Period Ratio 



Exterior Perturber 



- 3 




- 6 [ ; 

1.94 1.96 1.98 2 2.02 2.04 
Period Ratio 

Figure 4.7: This shows the perturber mass hmit from a more refined TTV analysis 
near the 2:1 mean- motion resonances as a function of the period ratio of the perturber 
and TrES-lb. The upper panel is for an interior perturber and the lower panel is for 
an exterior perturber. The solid curve is for a perturber eccentricity of 0.00, the 
dashed curve is for an eccentricity of 0.02, and the dotted curve is for an eccentricity 
of 0.05. The long-dashed line shows a mass of 5Mq, the dashed line is for 2M^, and 
the dotted line is for IM®. For a perturber with eccentricity 0.02, the lowest point 
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Figure 4.8: This shows the perturber mass hmit from a more refined TTV analysis 
near the 3:2 mean-motion resonances as a function of the period ratio of the perturber 
and TrES-lb. The legend for this plot is the same as in Figure WTfi 
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4.4 Discussion of TrES-1 Results 

Ultimately a combined analysis of all available data, including studies of the stability 
of candidate systems, will provide the most robust and sensitive determination of the 
parameters of each planetary system. None of the planetary systems that compose the 
limits in Figures [4.71 and Wl^ below lOM® are stable for more than 10^ orbits; though 
stable orbits with comparable masses, periods, and eccentricities exist. An overall 
stability analysis to accompany these TTV analyses was prohibitively expensive. I 
estimate that for Figure l^^ there were ~ 10^ potential systems that were analysed and 
~ 10^ systems were analysed for each curve in Figures 14.71 and 14.81 — requiring a total 
of ~ 10® hours of processor time. General stability analyses are often better suited 



to confirm or study a particular candidate system (e.g. ^, 



I, 3^) than to constrain 
candidate systems in the manner that we have done in this section. We did conduct 
a stability analysis for the system shown in Figure 14.11 and found it to be unstable 
within 10® orbits; no neighboring stable systems of comparable were found. 

It is unclear what fraction of probable planetary companions are excluded by our 
results. Aside from the fact that many of the known, multi-planetary systems are in 



mean-motion orbital resonance (e.g. GJ 876), recent works by |67f , [53, and [721 show 
that, under fairly general initial conditions, small planets are readily trapped in low- 
order, mean-motion resonances with a gas-giant as it interacts with a protoplanetary 
disk and migrates inward. These results may imply that resonant systems are common 
among multiple-planet systems. If this is true, then results like Figure 14.21 actually 
exclude a much larger fraction of the probable orbits than one might infer from the 



projection of the excluded regions onto the a/e plane. In addition, [4l| show that a 
large fraction of existing, terrestrial planets can survive the migration of a Jupiter- 
mass planet, though only a fraction of the survivors will be in resonance. 

The sensitivity of TTV to the mass of a perturbing planet renders it ideal for 
discovering and constraining the presence of additional planets in transiting systems 
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like TrES-1. These studies can help determine the ubiquity of multiple planet systems 
and resonant systems — including the distribution of mass in those systems. Moreover, 
TTV analyses of several systems can play a role in identifying the importance of 
various planet-formation mechanisms. For example, the presence of close-in terrestrial 
planets favors a sequential-accretion model of planet formation over a gravitational 
instability model |72l |. 

For the case of TrES-1, while a companion planet with a mass larger than the 
Earth is ruled out in a portion of the available parameter space, there remain large 
regions where additional planets could reside; this includes mutually inclined orbits 
which were not studied in this work (I believe that these results are valid provided 
that the angle of mutual inclination i satisfies cosz ~ 1). The 11 timings analysed in 
this work, excepting the questionable "12*^" point, deviate from a constant period by 
an amount that is difficult to interpret; having a of about 2 per degree of freedom. 
Additional observations, with higher precision, will allow us to more fully constrain 
the system and to interpret the existing transit timing variations. 



69 



Chapter 5 
THE HD 209458 SYSTEM 



We now turn our attention to the analysis of transit data from the HD209458 
system. The observations for HD209458 were reported by ^3 however, we re- 
reduced the raw data to obtain the transit times that are analyzed here. Many of the 
results of this chapter will be reported in the forthcoming P| . Since the data from this 
system were obtained using the Hubble Space Telescope, the precision of the timing 
measurements is significantly improved over those from the TrES-1 system. This 
makes these data capable of probing for yet smaller companions — well below the mass 
of the Earth in some regimes. Moreover, we will see that additional steps can be taken 
during observations to further improve the precision of the timing measurements. 

5.1 Data Reduction 

At least ten science observing programs with the Hubble Space Telescope have pro- 
posed to observe the transit of the planet in system HD209458. The five ultraviolet 
observations do not have enough photons to measure the transit lightcurve precisely, 
while five programs have made observations of the planetary transit in the optical and 
infrared, of which three are now publicly available: 8789, 9171 and 9447 (the other 
two will become available in the next year). 

Two programs. Brown 8789 and Charbonneau 9447, utilized the STIS spectro- 
graph as a CCD photometer, spreading the light from the star using a grism to 
capture as many photons as possible. This approach provides what may be the most 
precise relative photometry ever The third program utilized the Fine Guidance 
Sensors (FGS) to carry out rapid high-precision photometry jg^. These programs 
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Table 5.1: Transit epoch and number, instrument, wavelength range, cadence, and 
dates for HD209458 transit observations. 



Transit 


Orbs. 


Inst. 


Filt. 


wavelengths 


exp./read 


date 


(UT) 


1/0 


5 


STIS 


G750M 


581.3-638.2 nm 


60/20 sec 


2000 


Apr 25 


2/1 


5 


STIS 


G750M 


581.3-638.2 nm 


60/20 sec 


2000 


Apr 28-29 


3/3 


5 


STIS 


G750M 


581.3-638.2 nm 


60/20 sec 


2000 


May 5-6 


4/5 


5 


STIS 


G750M 


581.3-638.2 nm 


60/20 sec 


2000 


May 12-13 


5 /117 


3 


FGS 


F550W 


510-587.5 nm 


0.025/0 sec 


2001 


Jun 11 


6 /143 


2 


FGS 


F550W 


510-587.5 nm 


0.025/0 sec 


2001 


Sep 11 


7 /160 


2 


FGS 


F550W 


510-587.5 nm 


0.025/0 sec 


2001 


Nov 10 


8 /179 


2 


FGS 


F550W 


510-587.5 nm 


0.025/0 sec 


2002 


Jan 16 


9 /252 


2 


FGS 


F550W 


510-587.5 nm 


0.025/0 sec 


2002 


Sep 30 


10/313 


5 


STIS 


G430L 


290-570 nm 


22/20 sec 


2003 


May 3 


11/321 


5 


STIS 


G750L 


524-1027 nm 


19/20 sec 


2003 


May 31 


12/328 


5 


STIS 


G430L 


290-570 nm 


22/20 sec 


2003 


Jun 25 


13/331 


5 


STIS 


G750L 


524-1027 nm 


19/20 sec 


2003 


Jul 5-6 



give a total of thirteen transits observed with HST that are now available. A sum- 
mary of these observations is given in table 15.11 The remaining programs, with data 
that are not public, will provide another nine transits for future studies. 

The primary science objectives for these observations were to constrain the plan- 
etary parameters, to detect absorption features in the planet's atmosphere, and to 
search for satellites of the planet. Besides these studies, the data are useful for timing 
each transit to measure a precise ephemeris [t^ and for a TTV analysis. However, in 
order to obtain the transit times from each set of observations we needed to re-reduce 
the pipeline calibrated photometric data that is available in the archive. This process 
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involved several steps; following in part the reduction procedures described in 
and 



6^ . I summarize the most important steps here for each data set. For a more 
detailed description of the data reduction procedure see where the reduction is 
described in some detail. 

5.1.1 Outlier Removal and Data Binning 

For STIS we carried out the minimum reduction necessary to obtain precise photom- 
etry. We subtracted cosmic-rays using a 5 sigma rejection of the time series for each 
pixel within an HST orbit, we binned the photons within 8 pixels of the spectrum 
peak (17 pixels wide) to derive the total counts for each frame, and assigned the 
midpoint of each exposure as the time for each frame. We discarded the first frame 
from each orbit and for the first transit observed in program 8789 we only used the 
data for which the spectrum centroid was greater than 4 pixels from the edge of the 
CCD [see discussion in^^; in addition we did not utilize the first and last columns 
of the CCD which show larger errors. The 1-sigma errors on the flux are taken from 
the pipeline calibrated errors which are summed in quadrature. 

For FGS we took the sum of the four photometers, discarded data which deviated 
more than 3 sigma from a fourth-order polynomial fit to the time series for each orbit 
and binned the data to 80 second bins (10 second bins gave the same results) with the 
time stamp at the center of each bin. With these time series we then fit a model for 
the transiting planet and the flux sensitivity variations of Hubble, which I describe 
later. 

5.1.2 Photometric Error bars 

The errors computed from the pipeline include counting noise and, for the CCD data, 
read noise. We compared the scatter in the data to the pipeline error bars and found 
that the scatter was larger than the pipeline errors for the FGS and second set of 
STIS data. This could not be attributed solely to cosmic rays and we were unable to 
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identify the additional source of noise. Therefore, we estimate the photometric errors 
for each individual transit by analyzing the out-of-transit data. 

To account for the fluctuations in sensitivity of HST when determining the size 
of the error bars, we fit a model to the out-of-transit data for each transit with two 
components: i) a 5th order polynomial fit as a function of the orbital phase of HST 
and ii) a linear fit as a function of time. The residuals did not correlate with the HST 
orbital phase so we assigned the greater of the standard deviation of the residuals 
for the polynomial fit or the statistical error. The first five transits had statistical 
errors nearly equal to the residuals, while the error bars of the remaining transits 
were increased by 20-50%. These changes to the photometric error bars were used 
throughout the remainder of the analysis. 

5.1.3 Limb darkening 

If the lightcurves for each transit were only afi'ected by Poisson noise and were con- 
tinuously sampled, then assumptions about the shape of the transit lightcurve would 
be unnecessary: the mid-point of transit could be found by assuming the transit 
hghtcurve was symmetric. However, because the Earth occults the Hubble Space 
Telescope, the data are affected by sensitivity variations due to temperature changes 
of the telescope and by gaps in the observations. We must, therefore, choose a model 
for the transit lightcurve in order to identify the midpoint of the transit — the "tran- 
sit time" . The primary uncertainty in the modeling is the shape of the stellar limb 
darkening. We modeled the limb-darkening of the star with two approaches: i) stellar 
atmosphere predictions and ii) quadratic limb-darkening. By comparing the results 
from these two approaches we test the robustness of the transit times under changes 
to our assumptions about limb-darkening. 

Despite the differences between these these two approaches, the measured times 
were nearly identical within the errors and yielded a similar reduced x^- However, 
the inferred planetary radius and orbital inclination differed by more than the errors 
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for each model. 

5.1.4 HST flux sensitivity variations 

The thermal fluctuations of HST during its orbit lead to changes in the sensitivity 
of each instrument at the 0.5% level. Fortunately these variations appear to be 
reproducible during successive orbits (in addition to a smaller linear drift between 

n 

successive orbits), but they are not reproducible on longer timescales. Following [13[ 
and we phased the observations of each transit to the HST orbital period {ph)- 
We used two models for the variation of the flux sensitivity with orbital phase and 
two models for the "secular" variation with time over each transit, combined together 
giving a total of four models for the background variations. 

We modeled the sensitivity variations as a function of the orbital phase of HST as 
either a fifth-order polynomial, X]^=o '^^'i^''*^^''' "where < < 27r is the orbital phase 
or as a linear and two harmonic functions Co0 + Ci sin(0) + C2 cos(0) + C3 sin(20) + 
C4cos(20). The secular variation with time we modeled either as a linear function, 
c^ftranif) +£5^, whcrc ftran IS the transit lightcurve model normalized to unity outside 
transit or as a constant coefficient for each transit, Ci^5 ftran, i{t) where i = 0, 12 labels 
each transit. 

The resulting lightcurve model depends on sixteen non-linear parameters (the 3 
planetary orbital parameters and 13 transit times). For the models of the HST flux 
sensitivity variations there are either 91 or 113 linear parameters (6-9 for the HST 
sensitivity for each transit). Coincidentally, the physically interesting parameters 
are non-linear, while the uninteresting parameters are linear. We mapped the as a 
function of each non-linear parameter (while marginalizing over the other parameters) 
and found that the has a quadratic shape near the minimum — indicating a unique 
solution. 
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5.1.5 Timing Error Computation 

The best-fit transit timings for the STIS data have reduced chi-square in the range 
— 1.07 — 1.25 for the four different models of the flux sensitivity. However, the 
reduced chi-squares for the FGS data (which contain fewer data points) were as high 
as 2.6 (for the 5th transit) indicating that the models were a poor fit to the data. 
We examined the residuals of the FGS data and found that some residuals within the 
transits were correlated. Thus, our best fitting models were not a proper description 
of the data — presumably due to either errors in the limb darkening or in the modeling 
of the HST flux sensitivity variations. In effort to incorporate these systematic effects 
into our estimate of the photometric error bars we used a Monte Carlo bootstrap 
simulation of the errors where, in order to maintain the correlations, we shifted the 
residuals for each transit by a random number of points. We also reversed the residuals 
for each orbit and repeated the shift. This process was repeated 200 times; each time 
adding the shifted residuals to the best-fit model and then re-fitting with the same 
procedure that was applied to the original data as described in the previous section. 
Our expectation is that this procedure will give a better estimate of the errors that 
are due to inaccuracies in our assumed flux sensitivity models. We found that the 
resulting errors were larger by a factor of a few than if we simply used Poisson errors 
for the bootstrap or resampled the residuals randomly. 

Although this process gives an indication of the size of the systematic errors, we 
also assessed the systematic errors by comparing the results from the four different 
HST flux sensitivity variation models. The differences in between the four differ- 
ent models was small — all had reduced near unity, so it is difficult to favor one 
particular sensitivity variation model. However, the differences between the derived 
transit times for each model was larger than the size of the error bars for a given 
model. We beheve that this discrepancy is probably due to imperfections in all of the 
four flux sensitivity models, so we have taken the mean and standard deviations of all 
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eight hundred Monte-Carlo simulations as an estimate of the transit times and their 
errors. Finally, to check the robustness of the stellar atmosphere limb-darkening we 
have used the quadratic limb-darkening law with the polynomial function of phase 
and linear secular term to fit the transit times. We find that these times agree within 
the errors with the times derived using the stellar atmosphere model and the same 
flux sensitivity model, so we conclude that the limb-darkening model provides an 
accurate fit to the data. 

To these measured transit times we apply a correction for the motion of the Earth 
about the barycenter of the solar system (a heliocentric correction is insufficient as 
it differs from barycentric by up to 5 seconds). The motion of the Hubble Space 
Telescope about the Earth contributes a negligible timing error. The resulting transit 
times are reported in Table 15.21 We find from the Monte Carlo simulations that there 
is no significant correlation between the individual transit times so we report the 
standard deviation of each time without the full covariance matrix. 

5.2 Ephemeris 

The accurate and realistic estimation of errors on the transit times allows for an 
accurate estimate of the ephemeris of the transiting planet. We find values of: Tq = 
2451659.93677±0.00009 BJD (error is 8 seconds) and P = 3.5247484±0.0000004 days 
(error is 0.03 seconds). The period and zero-point of the ephemeris are anticorrelated 
with a correlation coefficient of -0.8; this is due to the fact that an increased zero 
point requires the period to shrink to go through the later data points. The anti- 
correlation decreases to -0.2 if the zeroth eclipse is chosen to be halfway between the 
first and last eclipses. Given the estimated statistical plus systematic errors, the chi- 
square is = 6.0 for eleven degrees of freedem (thirteen transits minus two fitting 
parameters). This indicates that the transit times are consistent with being uniform 
and indicates that there may not be a second planet detectable with the current HST 
data. However, as with the TrES-1 system these data can be used to place constraints 
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Tabic 5.2: Transit epoch, transit time, and timing uncertainty for HD 209458 transit 
observations. 



Transit 


Transit time 


Error 


Error 




(- 2450000 BJUj 


(days) 


(seconds) 


1 


1659.93678 


0.00015 


12 


2 


1663.46150 


0.00025 


21 


3 


1670.51102 


0.00013 


11 


4 


1677.56044 


0.00037 


31 


5 


2072.33284 


0.00049 


42 


6 


2163.97574 


0.00028 


24 


7 


2223.89685 


0.00058 


50 


8 


2290.86741 


0.00051 


44 


9 


2548.17311 


0.00029 


24 


10 


2763.18306 


0.00018 


15 


11 


2791.38087 


0.00013 


10 


12 


2816.05431 


0.00021 


18 


13 


2826.62867 


0.00017 


14 
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Figure 5.1: This shows the deviations from a constant period of the transit times for 
the HD209458 system. The horizontal axis is in HJD — 2450000 and the vertical axis 
is in days. 
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on the masses of other planets in the system. 

5.3 Constraints on Secondary Planets 

I conducted an analysis similar to that for the TrES-1 system for the HD209458 
system using the transit times given in table 15. 2[ However, because the reduced 
for these transit times is less than unity, I did not attempt to find the best-fitting 
two-planet system since any improvement in the overall reduced would be suspect. 
The constraints on the mass of a secondary planet in the HD209458 system are shown 
in Figure 15.21 for an interior perturber and in Figure 15.31 for an exterior perturber. 

As expected (in light of the results for the TrES-1 system) near the low-order 
mean-motion resonances the maximum allowed mass is comparable to or less than 
the mass of the Earth. A comparison to the radial velocity technique, similar to 
Figure is shown in Figure where we see that both the increase in the number 
of data (13 instead of 11) and the increased timing precision (and average of 25s 
instead of 108s) renders the TTV technique superior to the radial velocity technique 
over a larger range of period ratios. 

5.4 Statements About Space and Ground Based Observations 

A comparison of the results for the TrES-1 system and for the HD209458 system 
demonstrate the power that the TTV detection technique has to constrain compan- 
ions in transiting systems using either ground-based observations or space-based ob- 
servations, or both. It is somewhat surprising that ground-based observations taken 
with modest telescopes are able to probe for such small planetary companions. The 
average timing uncertainty in the data from the TrES-1 system was 108s compared to 
25s for the HD209458 system. This factor of four in precision resulted in a comparable 
tightening of the constraints on the presence of secondary planets in the HD209458 
system (see Figure I^TSj) . However, the most precise of the data for the TrES-1 sys- 
tem, taken with the 1.2m telescope at the Whipple Observatory, had an uncertainty 
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Figure 5.2: The upper limit on the mass of an interior perturbing planet in the 
HD209458 system. The regions bounded by the contours correspond to a companion 
mass that is less than 100 (dotted), 10 (dashed), and 1 (solid) earth masses. 
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Figure 5.3: The upper limit on the mass of an exterior perturbing planet in the 
HD209458 system. The regions bounded by the contours correspond to a companion 
mass that is less than 100 (dotted), 10 (dashed), and 1 (sohd) earth masses. 
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Figure 5.4: A comparison of the limits expressed in this paper (sohd curve) and the 
expected hmits that are achieved with 13 randomly selected radial velocity measure- 
ments that have a precision of 5 m s^^ (dotted curve) and 1 m s^^ (dashed curve) 
for a secondary planet with an eccentricity of 0.05. The radial velocity curves are 
given by equation ()4.2p using Ax^ = 9, N = 13, Qo = 0, and Q = 2. Near the 4:3, 
3:4, and the 1:2 mean-motion resonances the TTV analysis of these data can place 
tighter constraints (about a factor of 3) on the mass of putative secondary planets 
than can the RV technique with 1 m s~^ precision. 
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Figure 5.5: This image shows the mass hmit of a secondary planet obtained from a 
TTV analysis of the 11 transits for the TrES-1 system (upper curve) with the limit 
obtained from the 13 transits of the HD209458 system (lower curve). 

of 26s — almost identical to the average of the space-based observations. Moreover, 
the ground-based observations were not optimized for determining the time of the 
transit. This shows that ground based observations can play an important role in 
the collection of transit data for TTV analysis. I note that the STIS observations, 
with an average uncertainty of 17s, consistently gave better timing precision than the 
FGS data and is therefore of higher quality than the best ground-based observations 
of TrES-1. However, the potential for high quality observations from the ground to 
probe for terrestrial planets remains. 

While space-based observations will almost always provide higher precision tran- 
sit times than ground-based observations, the availability of l-2m telescopes is sig- 
nificantly higher than the availability of the few space-based observatories. Thus, 
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ground-based transit observations will likely provide follow-up observations of plan- 
etary systems that are discovered by the CoRoT and Kepler missions; and entirely 
ground-based transit programs such as XO certainly have the capacity to provide 
data of sufficent quality for interesting TTV analyses. 

Another item that became clear during the reduction of the HD209458 transit 
data is that despite its lofty position above the Earth's atmosphere, it is apparent 
that the Hubble Space Telescope is not ideally suited to timing transits. The main 
problem is that its low-Earth orbit causes occultation of the Sun which leads to 
thermal variations within the orbit and causes occultation of the target star which 
leads to gaps and repointing. The thermal variations cause sensitivity variations, 
while the gaps make a model for limb-darkening necessary to measure the transit 
times. If the data errors behaved according the Poisson statistics, the expected time 
variations for the STIS data would be of order 3-5 seconds, while the derived times 
are 2-3 times larger. The FGS data had a decreased sensitivity relative to the STIS 
data due to a lower count rate and the observation of each transit for only two or 
three orbits, an observing strategy which conserved on Hubble time but was unsuited 
to the precise measurement of transit times. Future precise measurements of transit 
times are possibly better carried out with satellites in Earth-trailing orbit, such as 
the Spitzer Space Telescope. Spitzer has the extra advantage that the infrared is 
less subject to limb-darkening; however, it will collect fewer photons because of its 
smaller aperture and because it observes past the peak of the photon count spectrum 
of G-type stars. 
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Chapter 6 
CONCLUSION 



The results presented in this dissertation favor tlie thesis tliat by analyzing the 
variations in transit times of a transiting extrasolar planet one can extract important 
scientific information. This information may include the orbital elements of unseen, 
perturbing planets — or, constraints on the presence of such planets, the absolute mass 
and radii of the bodies in the system, and relevant facts that bear on the formation 
and evolution of planetary systems. Thus far the first TTV analysis, that of the 
transits of the TrES-1 system, constitutes the first, published probe for sub earth- 
mass planets orbiting a distant, main-sequence star. While no discovery was made, 
the constraints that it places on the presence of secondary planets in the important, 
low-order orbital resonances remains an important result because if the implications 
that it has for theories of planet formation [72] • The same can be said for the similar 
analysis — with similar results — of the transits from the HD209458 system. 

However, two systems with similar results does not a trend make. Before I can 
make definitive statements concerning the results from TTV analyses I need to analyze 
many more systems; systems that should result from existing and planned transit 
searches (e.g. j^, ^, At this time, much development remains before TTV can 
be systematically apphed to the large amounts of transit data that are expected over 
the next decade. 

The research that remains has observational, theoretical, and analysis aspects. 
Some of these topics include: 1) optimizing ground and space-based observations for 
the timing of planetary transits and characterizing the effects of coverage gaps and 
timing noise on the observational side; 2) theoretical development to understand the 
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nature of j:l mean- motion resonances, the effects of mutually inclined orbits, the 
long-term stability of multiple planet systems, and the consequences of additional 
planetary companions beyond two; and 3) optimizing the existing analysis software 
so that it runs more efficiently and that it is can be applied more readily to multiple 
systems instead of the current situation of one system at a time. 

As new data from transiting systems becomes available TTV should prove an im- 
portant tool in addressing some of the fundamental questions concerning planetary 
systems in general and in making important discoveries in multiple planet systems. 
This is in part because of the sensitivity of the TTV scheme to the mass of a per- 
turbing planet — it being capable of probing for terrestrial planets with data taken 
from modest, ground-based telescopes. It is also because of TTV's unique ability to 
break existing degeneracies in the parameters involved in planetary systems such as 
the mass/radius degeneracy mentioned earlier. Also, TTV is particularly well-suited 
to probe for small objects in mean-motion resonances and to possibly determine the 
relative inclinations of planetary systems. 



86 



BIBLIOGRAPHY 



[1] E. Agol, J. Steffen, R. Sari, and W. Clarkson. Mon. Not. Royal Astr. Soc, 
359:567, 2005. 

[2] E. Agol and J. H. Steffen. in preparation. 

[3] M. D. Albrow et al. ApJ, 535:176, 2000. 

[4] R. Alonso, T. Brown, G. Torres, D. Latham, A. Sozzetti, G. Mandushev, J. Bel- 
monte, D. Charbonneau, H. Deeg, E. Dunham, F. O'Donovan, and R. Stefanik. 
ApJ, 613:L153, 2004. 

[5] J. W. Barnes and J. J. Fortney. Transit Detectabihty of Ring Systems Around 
Extrasolar Giant Planets, astro-ph/0409506, September 2004. 

[6] R. Barnes and T. Quinn. ApJ, 611:494, 2004. 

[7] D. P. Bennett and others. The Galactic Exoplanet Survey Telescope (GEST). In 
Future EUV/UV and Visible Space Astrophysics Missions and Instrumentation. 
Edited by J. Chris Blades, Oswald H. W. Siegmund. Proceedings of the SPIE, 
Volume 4854, pp. 141-157 (2003)., pages 141-157, February 2003. 

[8] T. Borkovits, B. Erdi, E. Forgacs-Dajka, and T. Kovacs. On the detectabihty 
of long period perturbations in close hierarchical triple stellar systems. A&A, 
398:1091, 2003. 

[9] W. J. Borucki et al. SPIE, 4854:129, 2003. 

[10] W. J. Borucki et al. The Kepler mission: a wide-field-of-view photometer de- 
signed to determine the frequency of Earth-size planets around solar-like stars. In 
Future EUV/UV and Visible Space Astrophysics Missions and Instrumentation. 
Edited by J. Chris Blades, Oswald H. W. Siegmund. Proceedings of the SPIE, 
Volume 4854, pp. 129-140 (2003)., pages 129-140, February 2003. 

[11] F. Bouchy, F. Pont, N. C. Santos, C. Melo, M. Mayor, D. Queloz, and S. Udry. 
A&A, 421:L13, 2004. 

[12] P. Bourde, Rouan D., and Leger A. A&A, 406:P203, 2003. 

[13] T. Brown, Charbonneau D., Gilliland R., Noyes R., and Burrows A. ApJ, 
552:699, 2001. 

[14] T. M. Brown. Expected Detection and False Alarm Rates for Transiting Jovian 
Planets. ApJL, 593:L125-L128, August 2003. 



87 



[15] P. Butler, S. S. Vogt, G. W. Marcy, D. A. Fischer, J. T. Wright, G. W. Henry, 
G. Laughhn, and J. Lissauer. A Neptune-Mass Planet Orbiting the Nearby M 
Dwarf GJ 436. astro-ph/0408587, August 2004. 

[16] R. P. Butler, G. W. Marcy, D. A. Fischer, T. M. Brown, A. R. Contos, S. G. 
Korzennik, P. Nisenson, and R. W. Noyes. Evidence for Multiple Companions 
to V Andromedae. ApJ, 526:916-927, December 1999. 

[17] D. Charbonneau. astro-ph/ 03 12252. 2003. 

[18] D. Charbonneau, T. M. Brown, D. W. Latham, and M. Mayor. ApJ, 529:145, 
2000. 

[19] D. Charbonneau et al. ApJ, 626:523, 2005. 

[20] A. Gumming. Detectability of Extrasolar Planets in Radial Velocity Surveys. 

astro-ph/04 084 70, August 2004. 

[21] H. J. Deeg, L. R. Doyle V. P. Kozhevnikov, J. E. Blue, E. L. Martin, and 
J. Schneider. A&A, 358:L5, 2000. 

[22] A. R. Dobrovolskis and W. J. Borucki. American Astronomical Society Division 

of Planetary Science Meeting, 28, 1208. 1996. 

[23] D. A. Fischer. G. W. Marcy. R. P. Butler, S. S. Vogt. S. Frink, and K. Apps. 
Planetary Companions to HD 12661, HD 92788, and HD 38529 and Variations 
in Keplerian Residuals of Extrasolar Planets. ApJ, 551:1107-1118, April 2001. 

[24] E. B. Ford and S. Trcmainc. PASP, 115, 2003. 

[25] A. Glindemann et al. Ap&SS, 286:35-44, 2003. 

[26] N. Gouda, T. Yano, Y. Yamada, Y. Kobayashi, T. Tsujimoto, and The Jasmine 
Working Group. In Proceedings of the GAIA Symposium, volume ESA-SP576, 
page 77, 2005. 

[27] A. Gould, B. S. Gaudi, and C. Han. Prospects for the Detection of Earth-Mass 
Planets. ArXiv Astrophysics e-prints. May 2004. 

[28] G. W. Henry, G. W. Marcy, R. P. Butler, and S. S. Vogt. ApJ, 529:L41, 2000. 

[29] M. Holman and N. Murray. Science, 307(5713) :1288, 2005. 

[30] K. Horne. Status aand Prospects of Planetary Transit Searches: Hot Jupiters 
Galore. In ASP Conf. Ser. 294-' Scientific Frontiers in Research on Eodrasolar 
Planets, page 361, 2003. 

[31] S. B. Howell, M. E. Everett, J. L. Tonry, A. Pickles, and C. Dain. Photometric 
Observations Using Orthogonal Transfer CCDs. PASP, 115:1340-1350, Novem- 
ber 2003. 



88 



[32] J. Ji, L. Liu, H. Kinoshita, and G. Li. ApJ, 631:1191, 2005. 

[33] M. Konacki, G. Torres, S. Jha, and D. D. Sassclov. An extrasolar planet that 
transits the disk of its parent star. Nature, 421:507-509, January 2003. 

[34] M. Konacki, G. Torres, D. D. Sasselov, G. Pietrzynski, A. Udalski, S. Jha, M. T. 
Ruiz, W. Gieren, and D. Minniti. A transiting extrasolar giant planet around 
the star OGLE-TR-113. astro-ph/0404541, April 2004. 

[35] M. G. Lattanzi et al. In Proceedings of the GAIA Symposium, volume ESA- 
SP576, page 251, 2005. 

[36] G. Laughhn, R. P. Butler, D. A. Fischer, G. W. Marcy, S. S. Vogt, and A. S. 
Wolf. The GJ 876 Planetary System - A Progress Report. astro-ph/040744i, 
July 2004. 

[37] G. Laughlin and J. E. Chambers. Short-Term Dynamical Interactions among 
Extrasolar Planets. ApJ, 551:L109, 2001. 

[38] R. Malhotra. in Planets Around Pulsars, ASP Conference Proceedings, J. A. 
Phillips, J. E. Thorsett, and S. R. Kulkarni, eds., page 89, 1993. 

[39] R. Malhotra. ApJ, 407:266, 1993. 

[40] K. Mandel and E. Agol. Analytic Lightcurves for Planetary Transit Searches. 
ApJ, 529:L45, 2002. 

[41] A. M. Mandell and S. Sigurdsson. ApJ, 599:L111, 2003. 

[42] G. W. Marcy, R. P. Butler, D. Fischer, S. S. Vogt, J. J. Lissauer, and E. J. 
Rivera. A Pair of Resonant Planets Orbiting GJ 876. ApJ, 556:296-301, July 
2001. 

[43] G. W. Marcy, R. P. Butler, D. A. Fischer, G. Laughhn, S. S. Vogt, G. W. Henry, 
and D. Pourbaix. A Planet at 5 AU around 55 Cancri. ApJ, 581:1375-1388, 
December 2002. 

[44] M. Mayor and D. Queloz. Nature, 378:355, 1995. 

[45] T. Mazeh et al. ApJ, 532:L55, 2000. 

[46] B. E. McArthur, M. Endl, W. D. Cochran, G. F. Benedict, D. A. Fischer, G. W. 
Marcy, R. P. Butler, D. Naef, M. Mayor, D. Queloz, S. Udry, and T. E. Harrison. 
Detection of a NEPTUNE-mass planet in the $\rho"{l}$ Cancri system using 
the Hobby-Eberly Telescope, astro-ph/0408585, August 2004. 

[47] P. R. McCuHough, J. E. Stys, J. A. Valenti, S. W. Fleming, K. A. Janes, and 
J. N. Heasley. PASP, 117,:783, 2005. 



89 



[48] J. Miralda-Escude. Orbital Perturbations of Transiting Planets: A Possible 
Method to Measure Stellar Quadrupoles and to Detect Earth-Mass Planets. ApJ, 
564:1019, 2002. 

[49] C. D. Murray and S. F. Dermott. Solar System Dynamics. (Cambridge: Cam- 
bridge Univ. Press), 1999. 

[50] M. Nauenberg. Determination of Masses and Other Properties of Extrasolar 
Planetary Systems with More than One Planet. ApJ, 568:369-376, March 2002. 

[51] M. Pan and R. Sari. A Generalization of the Lagrangian Points: Studies of 
Resonance for Highly Eccentric Orbits. AJ, 128:1418-1429, September 2004. 

[52] J. C. B. Papaloizou and E. Szuszkiewicz. Mon. Not. Royal Astr. Soc. Accepted, 
2005. 

[53] M. A. C. Perryman. Rep. Prog. Phys., 63:1209, 2000. 

[54] M. A. C. Perryman, K. S. dc Boer, G. Gilmorc, E. H0g, M. G. Lattanzi, L. Lin- 
degren, X. Luri, F. Mignard, O. Pace, and P. T. de Zeeuw. GAIA: Composition, 
formation and evolution of the Galaxy. A&A, 369:339-363, April 2001. 

[55] F. Pont, F. Bouchy, D. Queloz, N. C. Santos, C. Melo, M. Mayor, and S. Udry. 
A&A, 426:115, 2004. 

[56] S. H. Pravdo, S. B. Shaklan, T. Henry, and G. F. Benedict. ApJ, 617:1323, 2004. 

[57] W. H. Press, S. A. Tcukolsky, W. T. Vetterling, and B. P. Flannery. Numerical 
Recipes in C: the art of scientific computing. Cambridge University Press, second 
edition edition, 1992. 

[58] R. R. Rafikov. AJ, 125:922, 2003. 

[59] E. J. Rivera and J. J. Lissaucr. Dynamical Models of the Resonant Pair of 
Planets Orbiting the Star GJ 876. ApJ, 558:392-402, September 2001. 

[60] N. C. Santos, F. Bouchy, M. Mayor, F. Pepe, D. Queloz, S. Udry, C. Lovis, 
M. Bazot, W. Benz, J. . Bertaux, G. L. Curto, X. Delfosse, C. Mordasini, D. Naef, 
J. . Sivan, and S. Vauclair. The HARPS survey for southern extra-solar planets 
II. A 14 Earth-masses exoplanet around mu Arae. astro-ph/0408471, August 
2004. 

[61] P. Sartoretti and J. Schneider. On the detection of sateUites of extrasolar planets 
with the method of transits. A&As, 134:553-560, February 1999. 

[62] A. B. Schultz et al. Hst/fgs photometry of planetary transits of hd 209458. In 
AIP Conf. Proc. 713: The Search for Other Worlds, page 161, 2004. 

[63] S. Seager and G. Mallcn-Ornelas. ApJ, 585, 2003. 



90 



[64] S. Seagroves, J. Harkcr, G. Laughlin, J. Lacy, and T. Castellano. Detection 
of Intermediate-Period Transiting Planets with a Network of Small Telescopes: 
transitsearch.org. PASP, 115:1355-1362, December 2003. 

[65] A. Sozzetti, S. Casertano, R. A. Brown, and M. G. Lattanzi. PASP, 115, 2003. 

[66] J. H. Steffen and E. Agol. Mon. Not. Royal Astr. Soc, 364:L96, 2005. 

[67] E. Thommes. ApJ, 626:1033, 2005. 

[68] Y. Tsapras, K. Home, S. Kane, and C. Richard. MNRAS, 343:1131, 2003. 

[69] S. S. Vogt, G. W. Marcy, R. P. Butler, and K. Apps. Six New Planets from the 
Keck Precision Velocity Survey. ApJ, 536:902-914, June 2000. 

[70] J. Wisdom. The resonance overlap criterion and the onset of stochastic behavior 
in the restricted three-body problem. AJ, 85:1122-1133, August 1980. 

[71] R. A. Wittenmyer et al. ApJ, 632:1157, 2005. 

[72] J. L. Zhou, S. J. Aarseth, D. N. C. Lin, and M. Nagasawa. ApJ, 631:L85, 2005. 



91 



VITA 

Jason Hyrum Steffen was born in Fairfield, California on May 15, 1975 to David 
and Thail Steffen. He moved to Price, Utah when he was a few months old. After 
a brief excursion to Orangeville, Utah he moved to Fruit Heights, Utah where he 
attended grade through high school — graduating in 1993 with High Honors from Davis 
High School where he participated in almost everything that related to music and 
percussion. He then served a two-year mission for the Church of Jesus Christ of 
Latter-Day Saints in Brussels, Belgium. The day prior to his return to the United 
States, he had a brief, uneventful conversation with one Faith Kimball who, two 
years later, he met again. A more lengthy series of eventful conversations ensued and 
they were married in Salt Lake City, Utah on March 11, 2000. On May 16, 2005 
in Kirkland, Washington they were joined by James Hyrum Steffen, an invaluable 
addition to the family. 

Jason attended Weber State University in Ogden, Utah where he graduated 
Summa Cum Laude in Physics and Mathematics in June of 2000. For a year he 
worked as a Pricing Analyst for McLeod USA, an Embedded Software Engineer for 
L-3 Communications, and as a part-time faculty member of the Mathematics depart- 
ment at Salt Lake Community College. In 2001 Jason enrolled in the graduate school 
at the University of Washington where he studied Physics with Julianne Dalcanton, 
Paul Boynton, and finally Eric Agol. While there he occasionally taught Physics and 
Astronomy at Edmonds and Bellevue Community Colleges. His plan upon comple- 
tion of the requirements for his Doctor of Philosophy in Physics is to move to the 
Chicago area where he accepted the Brinson Postdoctoral Fellowship offered by the 
Particle Astrophysics Center at the Fermi National Accelerator Laboratory. 



